source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceerr2.F90 @ 8410

Last change on this file since 8410 was 8410, checked in by clem, 3 years ago

correction on the last commit

File size: 7.9 KB
Line 
1MODULE iceerr2
2   !!======================================================================
3   !!                     ***  MODULE  iceerr2  ***
4   !!   LIM-3 : Update of sea-ice global variables at the end of the time step
5   !!======================================================================
6   !! History :  3.0  !  2006-04  (M. Vancoppenolle) Original code
7   !!            3.5  !  2014-06  (C. Rousset)       Complete rewriting/cleaning
8   !!----------------------------------------------------------------------
9#if defined key_lim3
10   !!----------------------------------------------------------------------
11   !!   'key_lim3'                                      LIM3 sea-ice model
12   !!----------------------------------------------------------------------
13   !!    ice_err2   : computes update of sea-ice global variables from trend terms
14   !!----------------------------------------------------------------------
15   USE dom_oce
16   USE phycst          ! physical constants
17   USE ice
18   USE thd_ice         ! LIM thermodynamic sea-ice variables
19   USE limitd_th
20   USE limvar
21   USE limcons         ! conservation tests
22   USE limctl
23   !
24   USE lbclnk          ! lateral boundary condition - MPP exchanges
25   USE lib_mpp         ! MPP library
26   USE lib_fortran     ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
27   USE in_out_manager
28   USE timing          ! Timing
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC   ice_err2   ! routine called by ice_step
34
35   !! * Substitutions
36#  include "vectopt_loop_substitute.h90"
37   !!----------------------------------------------------------------------
38   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
39   !! $Id: iceerr2.F90 8378 2017-07-26 13:55:59Z clem $
40   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
41   !!----------------------------------------------------------------------
42CONTAINS
43
44   SUBROUTINE ice_err2( kt )
45      !!-------------------------------------------------------------------
46      !!               ***  ROUTINE ice_err2  ***
47      !!               
48      !! ** Purpose :  Computes update of sea-ice global variables at
49      !!               the end of the time step.
50      !!
51      !!---------------------------------------------------------------------
52      INTEGER, INTENT(in) ::   kt    ! number of iteration
53      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
54      REAL(wp) ::   zsal
55      REAL(wp) ::   zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 
56      !!-------------------------------------------------------------------
57      IF( nn_timing == 1 )  CALL timing_start('iceerr2')
58
59      IF( kt == nit000 .AND. lwp ) THEN
60         WRITE(numout,*)''
61         WRITE(numout,*)' ice_err2 '
62         WRITE(numout,*)' ~~~~~~~~~~~ '
63      ENDIF
64
65      ! conservation test
66      IF( ln_limdiachk ) CALL lim_cons_hsm(0, 'iceerr2', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b)
67
68      !----------------------------------------------------------------------
69      ! Constrain the thickness of the smallest category above himin
70      !----------------------------------------------------------------------
71      DO jj = 1, jpj 
72         DO ji = 1, jpi
73            rswitch = MAX( 0._wp , SIGN( 1._wp, a_i(ji,jj,1) - epsi20 ) )   !0 if no ice and 1 if yes
74            ht_i(ji,jj,1) = v_i (ji,jj,1) / MAX( a_i(ji,jj,1) , epsi20 ) * rswitch
75            IF( v_i(ji,jj,1) > 0._wp .AND. ht_i(ji,jj,1) < rn_himin ) THEN
76               a_i (ji,jj,1) = a_i (ji,jj,1) * ht_i(ji,jj,1) / rn_himin
77            ENDIF
78         END DO
79      END DO
80     
81      !-----------------------------------------------------
82      ! ice concentration should not exceed amax
83      !-----------------------------------------------------
84      at_i(:,:) = 0._wp
85      DO jl = 1, jpl
86         at_i(:,:) = a_i(:,:,jl) + at_i(:,:)
87      END DO
88
89      DO jl  = 1, jpl
90         DO jj = 1, jpj
91            DO ji = 1, jpi
92               IF( at_i(ji,jj) > rn_amax_2d(ji,jj) .AND. a_i(ji,jj,jl) > 0._wp ) THEN
93                  a_i (ji,jj,jl) = a_i (ji,jj,jl) * ( 1._wp - ( 1._wp - rn_amax_2d(ji,jj) / at_i(ji,jj) ) )
94               ENDIF
95            END DO
96         END DO
97      END DO
98
99      !---------------------
100      ! Ice salinity
101      !---------------------
102      IF (  nn_icesal == 2  ) THEN
103         DO jl = 1, jpl
104            DO jj = 1, jpj 
105               DO ji = 1, jpi
106                  zsal            = smv_i(ji,jj,jl)
107                  ! salinity stays in bounds
108                  rswitch         = 1._wp - MAX( 0._wp, SIGN( 1._wp, - v_i(ji,jj,jl) ) )
109                  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) )
110                  ! associated salt flux
111                  sfx_res(ji,jj) = sfx_res(ji,jj) - ( smv_i(ji,jj,jl) - zsal ) * rhoic * r1_rdtice
112               END DO
113            END DO
114         END DO
115      ENDIF
116
117      !----------------------------------------------------
118      ! Rebin categories with thickness out of bounds
119      !----------------------------------------------------
120      IF ( jpl > 1 ) CALL lim_itd_th_reb
121
122      !-----------------
123      ! zap small values
124      !-----------------
125      CALL lim_var_zapsmall
126
127      !------------------------------------------------------------------------------
128      ! Corrections to avoid wrong values                                        |
129      !------------------------------------------------------------------------------
130      ! Ice drift
131      !------------
132      DO jj = 2, jpjm1
133         DO ji = 2, jpim1
134            IF ( at_i(ji,jj) == 0._wp ) THEN ! what to do if there is no ice
135               IF ( at_i(ji+1,jj) == 0._wp ) u_ice(ji,jj)   = 0._wp ! right side
136               IF ( at_i(ji-1,jj) == 0._wp ) u_ice(ji-1,jj) = 0._wp ! left side
137               IF ( at_i(ji,jj+1) == 0._wp ) v_ice(ji,jj)   = 0._wp ! upper side
138               IF ( at_i(ji,jj-1) == 0._wp ) v_ice(ji,jj-1) = 0._wp ! bottom side
139            ENDIF
140         END DO
141      END DO
142      !lateral boundary conditions
143      CALL lbc_lnk_multi( u_ice(:,:), 'U', -1., v_ice(:,:), 'V', -1. )
144      !mask velocities
145      u_ice(:,:) = u_ice(:,:) * umask(:,:,1)
146      v_ice(:,:) = v_ice(:,:) * vmask(:,:,1)
147 
148      ! -------------------------------------------------
149      ! Diagnostics
150      ! -------------------------------------------------
151      DO jl  = 1, jpl
152         oa_i(:,:,jl) = oa_i(:,:,jl) + a_i(:,:,jl) * rdt_ice          ! ice natural aging
153         afx_thd(:,:) = afx_thd(:,:) + ( a_i(:,:,jl) - a_i_b(:,:,jl) ) * r1_rdtice
154      END DO
155      afx_tot = afx_thd + afx_dyn
156
157      DO jj = 1, jpj
158         DO ji = 1, jpi           
159            ! heat content variation (W.m-2)
160            diag_heat(ji,jj) = diag_heat(ji,jj) -  &
161               &               ( SUM( e_i(ji,jj,1:nlay_i,:) - e_i_b(ji,jj,1:nlay_i,:) ) +  & 
162               &                 SUM( e_s(ji,jj,1:nlay_s,:) - e_s_b(ji,jj,1:nlay_s,:) )    &
163               &               ) * r1_rdtice   
164            ! salt, volume
165            diag_smvi(ji,jj) = diag_smvi(ji,jj) + SUM( smv_i(ji,jj,:) - smv_i_b(ji,jj,:) ) * rhoic * r1_rdtice
166            diag_vice(ji,jj) = diag_vice(ji,jj) + SUM( v_i  (ji,jj,:) - v_i_b  (ji,jj,:) ) * rhoic * r1_rdtice
167            diag_vsnw(ji,jj) = diag_vsnw(ji,jj) + SUM( v_s  (ji,jj,:) - v_s_b  (ji,jj,:) ) * rhosn * r1_rdtice
168         END DO
169      END DO
170
171      ! conservation test
172      IF( ln_limdiachk ) CALL lim_cons_hsm(1, 'iceerr2', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b)
173
174      ! control prints
175      IF( ln_limctl )    CALL lim_prt( kt, iiceprt, jiceprt, 2, ' - Final state - ' )
176      IF( ln_ctl )       CALL lim_prt3D( 'iceerr2' )
177   
178      IF( nn_timing == 1 )  CALL timing_stop('iceerr2')
179
180   END SUBROUTINE ice_err2
181
182#else
183   !!----------------------------------------------------------------------
184   !!   Default option         Empty Module               No sea-ice model
185   !!----------------------------------------------------------------------
186CONTAINS
187   SUBROUTINE ice_err2     ! Empty routine
188   END SUBROUTINE ice_err2
189
190#endif
191
192END MODULE iceerr2
Note: See TracBrowser for help on using the repository browser.