source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icecor.F90 @ 8517

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

changes in style - part6 - one more round

File size: 10.3 KB
Line 
1MODULE icecor
2   !!======================================================================
3   !!                     ***  MODULE  icecor  ***
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_cor      : computes update of sea-ice global variables from trend terms
14   !!----------------------------------------------------------------------
15   USE dom_oce        ! ocean domain
16   USE phycst         ! physical constants
17   USE ice            ! sea-ice: variable
18   USE ice1D          ! sea-ice: thermodynamic sea-ice variables
19   USE iceitd         ! sea-ice: rebining
20   USE icevar         ! sea-ice: operations
21   USE icectl         ! sea-ice: control prints
22   !
23   USE in_out_manager ! I/O manager
24   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined) 
25   USE lbclnk         ! lateral boundary condition - MPP link
26   USE lib_mpp        ! MPP library
27   USE timing         ! Timing
28   USE iom            !
29
30   IMPLICIT NONE
31   PRIVATE
32
33   PUBLIC   ice_cor   ! called by icestp.F90
34
35   !! * Substitutions
36#  include "vectopt_loop_substitute.h90"
37   !!----------------------------------------------------------------------
38   !! NEMO/ICE 4.0 , NEMO Consortium (2017)
39   !! $Id: icecor.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_cor( kt, kn )
45      !!----------------------------------------------------------------------
46      !!               ***  ROUTINE ice_cor  ***
47      !!               
48      !! ** Purpose :   Computes corrections on sea-ice global variables at
49      !!              the end of the dynamics.
50      !!----------------------------------------------------------------------
51      INTEGER, INTENT(in) ::   kt    ! number of iteration
52      INTEGER, INTENT(in) ::   kn    ! 1 = after dyn ; 2 = after thermo
53      !
54      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices
55      REAL(wp) ::   zsal, zzc
56      REAL(wp), DIMENSION(jpi,jpj) ::   zafx   ! concentration trends diag
57      !!----------------------------------------------------------------------
58      ! controls
59      IF( nn_timing == 1 )   CALL timing_start('icecor')                                                             ! timing
60      IF( ln_icediachk   )   CALL ice_cons_hsm(0, 'icecor', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) ! conservation
61      !
62      IF( kt == nit000 .AND. lwp .AND. kn == 2 ) THEN
63         WRITE(numout,*)
64         WRITE(numout,*) 'ice_cor:  correct sea ice variables if out of bounds ' 
65         WRITE(numout,*) '~~~~~~~'
66      ENDIF
67      !
68      IF( kn == 2 ) THEN
69         !                          !-----------------------------------------------------
70         !                          !  thickness of the smallest category above himin    !
71         !                          !-----------------------------------------------------
72         WHERE( a_i(:,:,1) >= epsi20 )   ;   ht_i(:,:,1) = v_i (:,:,1) / a_i(:,:,1)
73         ELSEWHERE                       ;   ht_i(:,:,1) = 0._wp
74         END WHERE
75         WHERE( ht_i(:,:,1) < rn_himin )     a_i (:,:,1) = a_i (:,:,1) * ht_i(:,:,1) / rn_himin
76         !
77      ENDIF
78      !                             !-----------------------------------------------------
79      !                             !  ice concentration should not exceed amax          !
80      !                             !-----------------------------------------------------
81      at_i(:,:) = SUM( a_i(:,:,:), dim=3 )
82      DO jl  = 1, jpl
83         WHERE( at_i(:,:) > rn_amax_2d(:,:) )   a_i(:,:,jl) = a_i(:,:,jl) * rn_amax_2d(:,:) / at_i(:,:)
84      END DO
85   
86      !                             !-----------------------------------------------------
87      IF ( nn_icesal == 2 ) THEN    !  salinity must stay in bounds [Simin,Simax]        !
88      !                             !-----------------------------------------------------
89         zzc = rhoic * r1_rdtice
90         DO jl = 1, jpl
91            DO jj = 1, jpj 
92               DO ji = 1, jpi
93                  zsal = smv_i(ji,jj,jl)
94                  smv_i(ji,jj,jl) = MIN(  MAX( rn_simin*v_i(ji,jj,jl) , smv_i(ji,jj,jl) ) , rn_simax*v_i(ji,jj,jl)  )
95                  sfx_res(ji,jj) = sfx_res(ji,jj) - ( smv_i(ji,jj,jl) - zsal ) * zzc   ! associated salt flux
96               END DO
97            END DO
98         END DO
99      ENDIF
100
101      !                             !-----------------------------------------------------
102      !                             !  Rebin categories with thickness out of bounds     !
103      !                             !-----------------------------------------------------
104      IF ( jpl > 1 )   CALL ice_itd_reb( kt )
105
106      !                             !-----------------------------------------------------
107      CALL ice_var_zapsmall         !  Zap small values                                  !
108      !                             !-----------------------------------------------------
109
110      !                             !-----------------------------------------------------
111      IF( kn == 2 ) THEN            !  Ice drift case: Corrections to avoid wrong values !
112         DO jj = 2, jpjm1           !-----------------------------------------------------
113            DO ji = 2, jpim1
114               IF ( at_i(ji,jj) == 0._wp ) THEN    ! what to do if there is no ice
115                  IF ( at_i(ji+1,jj) == 0._wp )   u_ice(ji  ,jj) = 0._wp   ! right side
116                  IF ( at_i(ji-1,jj) == 0._wp )   u_ice(ji-1,jj) = 0._wp   ! left side
117                  IF ( at_i(ji,jj+1) == 0._wp )   v_ice(ji,jj  ) = 0._wp   ! upper side
118                  IF ( at_i(ji,jj-1) == 0._wp )   v_ice(ji,jj-1) = 0._wp   ! bottom side
119               ENDIF
120            END DO
121         END DO
122         CALL lbc_lnk_multi( u_ice, 'U', -1., v_ice, 'V', -1. )            ! lateral boundary conditions
123      ENDIF
124
125!!gm I guess the trends are only out on demand
126!!   So please, only do this is it exite an iom_use of on a these variables
127!!   furthermore, only allocate the diag_ arrays in this case
128!!   and do the iom_put here so that it is only a local allocation
129!!gm
130      !                             !-----------------------------------------------------
131      SELECT CASE( kn )             !  Diagnostics                                       !
132      !                             !-----------------------------------------------------
133      CASE( 1 )                        !--- dyn trend diagnostics
134         !
135!!gm   here I think the number of ice cat is too small to use a SUM instruction...
136         DO jj = 1, jpj
137            DO ji = 1, jpi           
138               !                 ! heat content variation (W.m-2)
139               diag_heat(ji,jj) = - (  SUM( e_i(ji,jj,1:nlay_i,:) - e_i_b(ji,jj,1:nlay_i,:) )    & 
140                  &                  + SUM( e_s(ji,jj,1:nlay_s,:) - e_s_b(ji,jj,1:nlay_s,:) )  ) * r1_rdtice
141               !                 ! salt, volume
142               diag_smvi(ji,jj) = SUM( smv_i(ji,jj,:) - smv_i_b(ji,jj,:) ) * rhoic * r1_rdtice
143               diag_vice(ji,jj) = SUM( v_i  (ji,jj,:) - v_i_b  (ji,jj,:) ) * rhoic * r1_rdtice
144               diag_vsnw(ji,jj) = SUM( v_s  (ji,jj,:) - v_s_b  (ji,jj,:) ) * rhosn * r1_rdtice
145            END DO
146         END DO
147         !                       ! concentration tendency (dynamics)
148         zafx   (:,:) = SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice 
149         afx_tot(:,:) = zafx(:,:)
150         IF( iom_use('afxdyn') )   CALL iom_put( 'afxdyn' , zafx(:,:) )
151         !
152      CASE( 2 )                        !--- thermo trend diagnostics & ice aging
153         !
154         oa_i(:,:,:) = oa_i(:,:,:) + a_i(:,:,:) * rdt_ice   ! ice natural aging incrementation
155         !
156!!gm   here I think the number of ice cat is too small to use a SUM instruction...
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) - (  SUM( e_i(ji,jj,1:nlay_i,:) - e_i_b(ji,jj,1:nlay_i,:) )    & 
161                  &                                   + SUM( e_s(ji,jj,1:nlay_s,:) - e_s_b(ji,jj,1:nlay_s,:) )  ) * r1_rdtice
162               !                 ! salt, volume
163               diag_smvi(ji,jj) = diag_smvi(ji,jj) + SUM( smv_i(ji,jj,:) - smv_i_b(ji,jj,:) ) * rhoic * r1_rdtice
164               diag_vice(ji,jj) = diag_vice(ji,jj) + SUM( v_i  (ji,jj,:) - v_i_b  (ji,jj,:) ) * rhoic * r1_rdtice
165               diag_vsnw(ji,jj) = diag_vsnw(ji,jj) + SUM( v_s  (ji,jj,:) - v_s_b  (ji,jj,:) ) * rhosn * r1_rdtice
166            END DO
167         END DO
168         !                       ! concentration tendency (total + thermo)
169         zafx   (:,:) = SUM( a_i(:,:,:) - a_i_b(:,:,:), dim=3 ) * r1_rdtice
170         afx_tot(:,:) = afx_tot(:,:) + zafx(:,:)
171         IF( iom_use('afxthd') )   CALL iom_put( 'afxthd' , zafx(:,:) )
172         IF( iom_use('afxtot') )   CALL iom_put( 'afxtot' , afx_tot(:,:) )
173         !
174      END SELECT
175      !
176      ! controls
177      IF( ln_icediachk   )   CALL ice_cons_hsm(1, 'icecor', rdiag_v, rdiag_s, rdiag_t, rdiag_fv, rdiag_fs, rdiag_ft) ! conservation
178      IF( ln_ctl         )   CALL ice_prt3D   ('icecor')                                                             ! prints
179      IF( ln_icectl .AND. kn == 2 )   CALL ice_prt( kt, iiceprt, jiceprt, 2, ' - Final state - ' )                   ! prints
180      IF( nn_timing == 1 )   CALL timing_stop ('icecor')                                                             ! timing
181      !
182   END SUBROUTINE ice_cor
183
184#else
185   !!----------------------------------------------------------------------
186   !!   Default option           Dummy module      NO LIM 3.0 sea-ice model
187   !!----------------------------------------------------------------------
188#endif
189
190   !!======================================================================
191END MODULE icecor
Note: See TracBrowser for help on using the repository browser.