source: NEMO/trunk/src/ICE/icectl.F90 @ 12377

Last change on this file since 12377 was 12377, checked in by acc, 10 months ago

The big one. Merging all 2019 developments from the option 1 branch back onto the trunk.

This changeset reproduces 2019/dev_r11943_MERGE_2019 on the trunk using a 2-URL merge
onto a working copy of the trunk. I.e.:

svn merge —ignore-ancestry \

svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/trunk \
svn+ssh://acc@forge.ipsl.jussieu.fr/ipsl/forge/projets/nemo/svn/NEMO/branches/2019/dev_r11943_MERGE_2019 ./

The —ignore-ancestry flag avoids problems that may otherwise arise from the fact that
the merge history been trunk and branch may have been applied in a different order but
care has been taken before this step to ensure that all applicable fixes and updates
are present in the merge branch.

The trunk state just before this step has been branched to releases/release-4.0-HEAD
and that branch has been immediately tagged as releases/release-4.0.2. Any fixes
or additions in response to tickets on 4.0, 4.0.1 or 4.0.2 should be done on
releases/release-4.0-HEAD. From now on future 'point' releases (e.g. 4.0.2) will
remain unchanged with periodic releases as needs demand. Note release-4.0-HEAD is a
transitional naming convention. Future full releases, say 4.2, will have a release-4.2
branch which fulfills this role and the first point release (e.g. 4.2.0) will be made
immediately following the release branch creation.

2020 developments can be started from any trunk revision later than this one.

  • Property svn:keywords set to Id
File size: 39.3 KB
Line 
1MODULE icectl
2   !!======================================================================
3   !!                     ***  MODULE  icectl  ***
4   !!   sea-ice : controls and prints
5   !!======================================================================
6   !! History :  3.5  !  2015-01  (M. Vancoppenolle) Original code
7   !!            3.7  !  2016-10  (C. Rousset)       Add routine ice_prt3D
8   !!            4.0  !  2018     (many people)      SI3 [aka Sea Ice cube]
9   !!----------------------------------------------------------------------
10#if defined key_si3
11   !!----------------------------------------------------------------------
12   !!   'key_si3'                                       SI3 sea-ice model
13   !!----------------------------------------------------------------------
14   !!    ice_cons_hsm     : conservation tests on heat, salt and mass during a  time step (global)
15   !!    ice_cons_final   : conservation tests on heat, salt and mass at end of time step (global)
16   !!    ice_cons2D       : conservation tests on heat, salt and mass at each gridcell
17   !!    ice_ctl          : control prints in case of crash
18   !!    ice_prt          : control prints at a given grid point
19   !!    ice_prt3D        : control prints of ice arrays
20   !!----------------------------------------------------------------------
21   USE phycst         ! physical constants
22   USE oce            ! ocean dynamics and tracers
23   USE dom_oce        ! ocean space and time domain
24   USE ice            ! sea-ice: variables
25   USE ice1D          ! sea-ice: thermodynamics variables
26   USE sbc_oce        ! Surface boundary condition: ocean fields
27   USE sbc_ice        ! Surface boundary condition: ice   fields
28   !
29   USE in_out_manager ! I/O manager
30   USE iom            ! I/O manager library
31   USE lib_mpp        ! MPP library
32   USE lib_fortran    ! fortran utilities (glob_sum + no signed zero)
33   USE timing         ! Timing
34   USE prtctl         ! Print control
35
36   IMPLICIT NONE
37   PRIVATE
38
39   PUBLIC   ice_cons_hsm
40   PUBLIC   ice_cons_final
41   PUBLIC   ice_cons2D
42   PUBLIC   ice_ctl
43   PUBLIC   ice_prt
44   PUBLIC   ice_prt3D
45
46   ! thresold rates for conservation
47   !    these values are changed by the namelist parameter rn_icechk, so that threshold = zchk * rn_icechk
48   REAL(wp), PARAMETER ::   zchk_m   = 2.5e-7   ! kg/m2/s <=> 1e-6 m of ice per hour spuriously gained/lost
49   REAL(wp), PARAMETER ::   zchk_s   = 2.5e-6   ! g/m2/s  <=> 1e-6 m of ice per hour spuriously gained/lost (considering s=10g/kg)
50   REAL(wp), PARAMETER ::   zchk_t   = 7.5e-2   ! W/m2    <=> 1e-6 m of ice per hour spuriously gained/lost (considering Lf=3e5J/kg)
51   
52   !! * Substitutions
53#  include "do_loop_substitute.h90"
54   !!----------------------------------------------------------------------
55   !! NEMO/ICE 4.0 , NEMO Consortium (2018)
56   !! $Id$
57   !! Software governed by the CeCILL license (see ./LICENSE)
58   !!----------------------------------------------------------------------
59CONTAINS
60
61   SUBROUTINE ice_cons_hsm( icount, cd_routine, pdiag_v, pdiag_s, pdiag_t, pdiag_fv, pdiag_fs, pdiag_ft )
62      !!-------------------------------------------------------------------
63      !!                       ***  ROUTINE ice_cons_hsm ***
64      !!
65      !! ** Purpose : Test the conservation of heat, salt and mass for each ice routine
66      !!                     + test if ice concentration and volume are > 0
67      !!
68      !! ** Method  : This is an online diagnostics which can be activated with ln_icediachk=true
69      !!              It prints in ocean.output if there is a violation of conservation at each time-step
70      !!              The thresholds (zchk_m, zchk_s, zchk_t) determine violations
71      !!              For salt and heat thresholds, ice is considered to have a salinity of 10
72      !!              and a heat content of 3e5 J/kg (=latent heat of fusion)
73      !!-------------------------------------------------------------------
74      INTEGER         , INTENT(in)    ::   icount        ! called at: =0 the begining of the routine, =1  the end
75      CHARACTER(len=*), INTENT(in)    ::   cd_routine    ! name of the routine
76      REAL(wp)        , INTENT(inout) ::   pdiag_v, pdiag_s, pdiag_t, pdiag_fv, pdiag_fs, pdiag_ft
77      !!
78      REAL(wp) ::   zdiag_mass, zdiag_salt, zdiag_heat, &
79         &          zdiag_vmin, zdiag_amin, zdiag_amax, zdiag_eimin, zdiag_esmin, zdiag_smin
80      REAL(wp) ::   zvtrp, zetrp
81      REAL(wp) ::   zarea
82      !!-------------------------------------------------------------------
83      !
84      IF( icount == 0 ) THEN
85
86         pdiag_v = glob_sum( 'icectl',   SUM( v_i * rhoi + v_s * rhos, dim=3 ) * e1e2t )
87         pdiag_s = glob_sum( 'icectl',   SUM( sv_i * rhoi            , dim=3 ) * e1e2t )
88         pdiag_t = glob_sum( 'icectl', ( SUM( SUM( e_i, dim=4 ), dim=3 ) + SUM( SUM( e_s, dim=4 ), dim=3 ) ) * e1e2t )
89
90         ! mass flux
91         pdiag_fv = glob_sum( 'icectl',  &
92            &                         ( wfx_bog + wfx_bom + wfx_sum + wfx_sni + wfx_opw + wfx_res + wfx_dyn + wfx_lam + wfx_pnd + &
93            &                           wfx_snw_sni + wfx_snw_sum + wfx_snw_dyn + wfx_snw_sub + wfx_ice_sub + wfx_spr ) * e1e2t )
94         ! salt flux
95         pdiag_fs = glob_sum( 'icectl',  &
96            &                         ( sfx_bri + sfx_bog + sfx_bom + sfx_sum + sfx_sni + &
97            &                           sfx_opw + sfx_res + sfx_dyn + sfx_sub + sfx_lam ) * e1e2t )
98         ! heat flux
99         pdiag_ft = glob_sum( 'icectl',  &
100            &                         (   hfx_sum + hfx_bom + hfx_bog + hfx_dif + hfx_opw + hfx_snw  &
101            &                           - hfx_thd - hfx_dyn - hfx_res - hfx_sub - hfx_spr ) * e1e2t )
102
103      ELSEIF( icount == 1 ) THEN
104
105         ! -- mass diag -- !
106         zdiag_mass = ( glob_sum( 'icectl', SUM( v_i * rhoi + v_s * rhos, dim=3 ) * e1e2t ) - pdiag_v ) * r1_rdtice       &
107            &         + glob_sum( 'icectl', ( wfx_bog + wfx_bom + wfx_sum + wfx_sni + wfx_opw + wfx_res + wfx_dyn +       &
108            &                                 wfx_lam + wfx_pnd + wfx_snw_sni + wfx_snw_sum + wfx_snw_dyn + wfx_snw_sub + &
109            &                                 wfx_ice_sub + wfx_spr ) * e1e2t )                                           &
110            &         - pdiag_fv
111         !
112         ! -- salt diag -- !
113         zdiag_salt = ( glob_sum( 'icectl', SUM( sv_i * rhoi , dim=3 ) * e1e2t ) - pdiag_s ) * r1_rdtice  &
114            &         + glob_sum( 'icectl', ( sfx_bri + sfx_bog + sfx_bom + sfx_sum + sfx_sni +           &
115            &                                 sfx_opw + sfx_res + sfx_dyn + sfx_sub + sfx_lam ) * e1e2t ) &
116            &         - pdiag_fs
117         !
118         ! -- heat diag -- !
119         zdiag_heat = ( glob_sum( 'icectl', ( SUM(SUM(e_i, dim=4), dim=3) + SUM(SUM(e_s, dim=4), dim=3) ) * e1e2t ) - pdiag_t &
120            &         ) * r1_rdtice                                                                                           &
121            &         + glob_sum( 'icectl', (  hfx_sum + hfx_bom + hfx_bog + hfx_dif + hfx_opw + hfx_snw                      &
122            &                                - hfx_thd - hfx_dyn - hfx_res - hfx_sub - hfx_spr ) * e1e2t )                    &
123            &         - pdiag_ft
124
125         ! -- min/max diag -- !
126         zdiag_amax  = glob_max( 'icectl', SUM( a_i, dim=3 ) )
127         zdiag_vmin  = glob_min( 'icectl', v_i )
128         zdiag_amin  = glob_min( 'icectl', a_i )
129         zdiag_smin  = glob_min( 'icectl', sv_i )
130         zdiag_eimin = glob_min( 'icectl', SUM( e_i, dim=3 ) )
131         zdiag_esmin = glob_min( 'icectl', SUM( e_s, dim=3 ) )
132
133         ! -- advection scheme is conservative? -- !
134         zvtrp = glob_sum( 'icectl', ( diag_trp_vi * rhoi + diag_trp_vs * rhos ) * e1e2t ) ! must be close to 0 (only for Prather)
135         zetrp = glob_sum( 'icectl', ( diag_trp_ei        + diag_trp_es        ) * e1e2t ) ! must be close to 0 (only for Prather)
136
137         ! ice area (+epsi10 to set a threshold > 0 when there is no ice)
138         zarea = glob_sum( 'icectl', SUM( a_i + epsi10, dim=3 ) * e1e2t )
139
140         IF( lwp ) THEN
141            ! check conservation issues
142            IF( ABS(zdiag_mass) > zchk_m * rn_icechk_glo * zarea ) &
143               &                   WRITE(numout,*)   cd_routine,' : violation mass cons. [kg] = ',zdiag_mass * rdt_ice
144            IF( ABS(zdiag_salt) > zchk_s * rn_icechk_glo * zarea ) &
145               &                   WRITE(numout,*)   cd_routine,' : violation salt cons. [g]  = ',zdiag_salt * rdt_ice
146            IF( ABS(zdiag_heat) > zchk_t * rn_icechk_glo * zarea ) &
147               &                   WRITE(numout,*)   cd_routine,' : violation heat cons. [J]  = ',zdiag_heat * rdt_ice
148            ! check negative values
149            IF( zdiag_vmin  < 0. ) WRITE(numout,*)   cd_routine,' : violation v_i < 0         = ',zdiag_vmin
150            IF( zdiag_amin  < 0. ) WRITE(numout,*)   cd_routine,' : violation a_i < 0         = ',zdiag_amin
151            IF( zdiag_smin  < 0. ) WRITE(numout,*)   cd_routine,' : violation s_i < 0         = ',zdiag_smin
152            IF( zdiag_eimin < 0. ) WRITE(numout,*)   cd_routine,' : violation e_i < 0         = ',zdiag_eimin
153            IF( zdiag_esmin < 0. ) WRITE(numout,*)   cd_routine,' : violation e_s < 0         = ',zdiag_esmin
154            ! check maximum ice concentration
155            IF( zdiag_amax > MAX(rn_amax_n,rn_amax_s)+epsi10 .AND. cd_routine /= 'icedyn_adv' .AND. cd_routine /= 'icedyn_rdgrft' ) &
156               &                   WRITE(numout,*)   cd_routine,' : violation a_i > amax      = ',zdiag_amax
157            ! check if advection scheme is conservative
158            !    only check for Prather because Ultimate-Macho uses corrective fluxes (wfx etc)
159            !    so the formulation for conservation is different (and not coded)
160            !    it does not mean UM is not conservative (it is checked with above prints) => update (09/2019): same for Prather now
161            !IF( ln_adv_Pra .AND. ABS(zvtrp) > zchk_m * rn_icechk_glo * zarea .AND. cd_routine == 'icedyn_adv' ) &
162            !   &                   WRITE(numout,*)   cd_routine,' : violation adv scheme [kg] = ',zvtrp * rdt_ice
163         ENDIF
164         !
165      ENDIF
166
167   END SUBROUTINE ice_cons_hsm
168
169   SUBROUTINE ice_cons_final( cd_routine )
170      !!-------------------------------------------------------------------
171      !!                     ***  ROUTINE ice_cons_final ***
172      !!
173      !! ** Purpose : Test the conservation of heat, salt and mass at the end of each ice time-step
174      !!
175      !! ** Method  : This is an online diagnostics which can be activated with ln_icediachk=true
176      !!              It prints in ocean.output if there is a violation of conservation at each time-step
177      !!              The thresholds (zchk_m, zchk_s, zchk_t) determine the violations
178      !!              For salt and heat thresholds, ice is considered to have a salinity of 10
179      !!              and a heat content of 3e5 J/kg (=latent heat of fusion)
180      !!-------------------------------------------------------------------
181      CHARACTER(len=*), INTENT(in) ::   cd_routine    ! name of the routine
182      REAL(wp) ::   zdiag_mass, zdiag_salt, zdiag_heat
183      REAL(wp) ::   zarea
184      !!-------------------------------------------------------------------
185
186      ! water flux
187      ! -- mass diag -- !
188      zdiag_mass = glob_sum( 'icectl', ( wfx_ice + wfx_snw + wfx_spr + wfx_sub + diag_vice + diag_vsnw ) * e1e2t )
189
190      ! -- salt diag -- !
191      zdiag_salt = glob_sum( 'icectl', ( sfx + diag_sice ) * e1e2t )
192
193      ! -- heat diag -- !
194      ! clem: not the good formulation
195      !!zdiag_heat  = glob_sum( 'icectl', ( qt_oce_ai - qt_atm_oi + diag_heat + hfx_thd + hfx_dyn + hfx_res + hfx_sub + hfx_spr  &
196      !!   &                              ) * e1e2t )
197
198      ! ice area (+epsi10 to set a threshold > 0 when there is no ice)
199      zarea = glob_sum( 'icectl', SUM( a_i + epsi10, dim=3 ) * e1e2t )
200
201      IF( lwp ) THEN
202         IF( ABS(zdiag_mass) > zchk_m * rn_icechk_glo * zarea ) &
203            &                   WRITE(numout,*) cd_routine,' : violation mass cons. [kg] = ',zdiag_mass * rdt_ice
204         IF( ABS(zdiag_salt) > zchk_s * rn_icechk_glo * zarea ) &
205            &                   WRITE(numout,*) cd_routine,' : violation salt cons. [g]  = ',zdiag_salt * rdt_ice
206         !!IF( ABS(zdiag_heat) > zchk_t * rn_icechk_glo * zarea ) WRITE(numout,*) cd_routine,' : violation heat cons. [J]  = ',zdiag_heat * rdt_ice
207      ENDIF
208      !
209   END SUBROUTINE ice_cons_final
210
211   SUBROUTINE ice_cons2D( icount, cd_routine, pdiag_v, pdiag_s, pdiag_t, pdiag_fv, pdiag_fs, pdiag_ft )
212      !!-------------------------------------------------------------------
213      !!                       ***  ROUTINE ice_cons2D ***
214      !!
215      !! ** Purpose : Test the conservation of heat, salt and mass for each ice routine
216      !!                     + test if ice concentration and volume are > 0
217      !!
218      !! ** Method  : This is an online diagnostics which can be activated with ln_icediachk=true
219      !!              It stops the code if there is a violation of conservation at any gridcell
220      !!-------------------------------------------------------------------
221      INTEGER         , INTENT(in) ::   icount        ! called at: =0 the begining of the routine, =1  the end
222      CHARACTER(len=*), INTENT(in) ::   cd_routine    ! name of the routine
223      REAL(wp)        , DIMENSION(jpi,jpj), INTENT(inout) ::   pdiag_v, pdiag_s, pdiag_t, pdiag_fv, pdiag_fs, pdiag_ft
224      !!
225      REAL(wp), DIMENSION(jpi,jpj) ::   zdiag_mass, zdiag_salt, zdiag_heat, &
226         &                              zdiag_amin, zdiag_vmin, zdiag_smin, zdiag_emin !!, zdiag_amax 
227      INTEGER ::   jl, jk
228      LOGICAL ::   ll_stop_m = .FALSE.
229      LOGICAL ::   ll_stop_s = .FALSE.
230      LOGICAL ::   ll_stop_t = .FALSE.
231      CHARACTER(len=120) ::   clnam   ! filename for the output
232      !!-------------------------------------------------------------------
233      !
234      IF( icount == 0 ) THEN
235
236         pdiag_v = SUM( v_i  * rhoi + v_s * rhos, dim=3 )
237         pdiag_s = SUM( sv_i * rhoi             , dim=3 )
238         pdiag_t = SUM( SUM( e_i, dim=4 ), dim=3 ) + SUM( SUM( e_s, dim=4 ), dim=3 )
239
240         ! mass flux
241         pdiag_fv = wfx_bog + wfx_bom + wfx_sum + wfx_sni + wfx_opw + wfx_res + wfx_dyn + wfx_lam + wfx_pnd  +  &
242            &       wfx_snw_sni + wfx_snw_sum + wfx_snw_dyn + wfx_snw_sub + wfx_ice_sub + wfx_spr
243         ! salt flux
244         pdiag_fs = sfx_bri + sfx_bog + sfx_bom + sfx_sum + sfx_sni + sfx_opw + sfx_res + sfx_dyn + sfx_sub + sfx_lam 
245         ! heat flux
246         pdiag_ft =   hfx_sum + hfx_bom + hfx_bog + hfx_dif + hfx_opw + hfx_snw  & 
247            &       - hfx_thd - hfx_dyn - hfx_res - hfx_sub - hfx_spr
248
249      ELSEIF( icount == 1 ) THEN
250
251         ! -- mass diag -- !
252         zdiag_mass =   ( SUM( v_i * rhoi + v_s * rhos, dim=3 ) - pdiag_v ) * r1_rdtice                             &
253            &         + ( wfx_bog + wfx_bom + wfx_sum + wfx_sni + wfx_opw + wfx_res + wfx_dyn + wfx_lam + wfx_pnd + &
254            &             wfx_snw_sni + wfx_snw_sum + wfx_snw_dyn + wfx_snw_sub + wfx_ice_sub + wfx_spr )           &
255            &         - pdiag_fv
256         IF( MAXVAL( ABS(zdiag_mass) ) > zchk_m * rn_icechk_cel )   ll_stop_m = .TRUE.
257         !
258         ! -- salt diag -- !
259         zdiag_salt =   ( SUM( sv_i * rhoi , dim=3 ) - pdiag_s ) * r1_rdtice                                                  &
260            &         + ( sfx_bri + sfx_bog + sfx_bom + sfx_sum + sfx_sni + sfx_opw + sfx_res + sfx_dyn + sfx_sub + sfx_lam ) &
261            &         - pdiag_fs
262         IF( MAXVAL( ABS(zdiag_salt) ) > zchk_s * rn_icechk_cel )   ll_stop_s = .TRUE.
263         !
264         ! -- heat diag -- !
265         zdiag_heat =   ( SUM( SUM( e_i, dim=4 ), dim=3 ) + SUM( SUM( e_s, dim=4 ), dim=3 ) - pdiag_t ) * r1_rdtice &
266            &         + (  hfx_sum + hfx_bom + hfx_bog + hfx_dif + hfx_opw + hfx_snw                                & 
267            &            - hfx_thd - hfx_dyn - hfx_res - hfx_sub - hfx_spr )                                        &
268            &         - pdiag_ft
269         IF( MAXVAL( ABS(zdiag_heat) ) > zchk_t * rn_icechk_cel )   ll_stop_t = .TRUE.
270         !
271         ! -- other diags -- !
272         ! a_i < 0
273         zdiag_amin(:,:) = 0._wp
274         DO jl = 1, jpl
275            WHERE( a_i(:,:,jl) < 0._wp )   zdiag_amin(:,:) = 1._wp
276         ENDDO
277         ! v_i < 0
278         zdiag_vmin(:,:) = 0._wp
279         DO jl = 1, jpl
280            WHERE( v_i(:,:,jl) < 0._wp )   zdiag_vmin(:,:) = 1._wp
281         ENDDO
282         ! s_i < 0
283         zdiag_smin(:,:) = 0._wp
284         DO jl = 1, jpl
285            WHERE( s_i(:,:,jl) < 0._wp )   zdiag_smin(:,:) = 1._wp
286         ENDDO
287         ! e_i < 0
288         zdiag_emin(:,:) = 0._wp
289         DO jl = 1, jpl
290            DO jk = 1, nlay_i
291               WHERE( e_i(:,:,jk,jl) < 0._wp )   zdiag_emin(:,:) = 1._wp
292            ENDDO
293         ENDDO
294         ! a_i > amax
295         !WHERE( SUM( a_i, dim=3 ) > ( MAX( rn_amax_n, rn_amax_s ) + epsi10 )   ;   zdiag_amax(:,:) = 1._wp
296         !ELSEWHERE                                                             ;   zdiag_amax(:,:) = 0._wp
297         !END WHERE
298
299         IF( ll_stop_m .OR. ll_stop_s .OR. ll_stop_t ) THEN
300            clnam = 'diag_ice_conservation_'//cd_routine
301            CALL ice_cons_wri( clnam, zdiag_mass, zdiag_salt, zdiag_heat, zdiag_amin, zdiag_vmin, zdiag_smin, zdiag_emin )
302         ENDIF
303
304         IF( ll_stop_m )   CALL ctl_stop( 'STOP', cd_routine//': ice mass conservation issue' )
305         IF( ll_stop_s )   CALL ctl_stop( 'STOP', cd_routine//': ice salt conservation issue' )
306         IF( ll_stop_t )   CALL ctl_stop( 'STOP', cd_routine//': ice heat conservation issue' )
307         
308      ENDIF
309
310   END SUBROUTINE ice_cons2D
311
312   SUBROUTINE ice_cons_wri( cdfile_name, pdiag_mass, pdiag_salt, pdiag_heat, pdiag_amin, pdiag_vmin, pdiag_smin, pdiag_emin )
313      !!---------------------------------------------------------------------
314      !!                 ***  ROUTINE ice_cons_wri  ***
315      !!       
316      !! ** Purpose :   create a NetCDF file named cdfile_name which contains
317      !!                the instantaneous fields when conservation issue occurs
318      !!
319      !! ** Method  :   NetCDF files using ioipsl
320      !!----------------------------------------------------------------------
321      CHARACTER(len=*), INTENT( in ) ::   cdfile_name      ! name of the file created
322      REAL(wp), DIMENSION(:,:), INTENT( in ) ::   pdiag_mass, pdiag_salt, pdiag_heat, &
323         &                                        pdiag_amin, pdiag_vmin, pdiag_smin, pdiag_emin !!, pdiag_amax 
324      !!
325      INTEGER ::   inum
326      !!----------------------------------------------------------------------
327      !
328      IF(lwp) WRITE(numout,*)
329      IF(lwp) WRITE(numout,*) 'ice_cons_wri : single instantaneous ice state'
330      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~  named :', cdfile_name, '...nc'
331      IF(lwp) WRITE(numout,*)               
332
333      CALL iom_open( TRIM(cdfile_name), inum, ldwrt = .TRUE., kdlev = jpl )
334     
335      CALL iom_rstput( 0, 0, inum, 'cons_mass', pdiag_mass(:,:) , ktype = jp_r8 )    ! ice mass spurious lost/gain
336      CALL iom_rstput( 0, 0, inum, 'cons_salt', pdiag_salt(:,:) , ktype = jp_r8 )    ! ice salt spurious lost/gain
337      CALL iom_rstput( 0, 0, inum, 'cons_heat', pdiag_heat(:,:) , ktype = jp_r8 )    ! ice heat spurious lost/gain
338      ! other diags
339      CALL iom_rstput( 0, 0, inum, 'aneg_count', pdiag_amin(:,:) , ktype = jp_r8 )    !
340      CALL iom_rstput( 0, 0, inum, 'vneg_count', pdiag_vmin(:,:) , ktype = jp_r8 )    !
341      CALL iom_rstput( 0, 0, inum, 'sneg_count', pdiag_smin(:,:) , ktype = jp_r8 )    !
342      CALL iom_rstput( 0, 0, inum, 'eneg_count', pdiag_emin(:,:) , ktype = jp_r8 )    !
343     
344      CALL iom_close( inum )
345
346   END SUBROUTINE ice_cons_wri
347   
348   SUBROUTINE ice_ctl( kt )
349      !!-------------------------------------------------------------------
350      !!                   ***  ROUTINE ice_ctl ***
351      !!                 
352      !! ** Purpose :   Alerts in case of model crash
353      !!-------------------------------------------------------------------
354      INTEGER, INTENT(in) ::   kt      ! ocean time step
355      INTEGER  ::   ji, jj, jk,  jl   ! dummy loop indices
356      INTEGER  ::   inb_altests       ! number of alert tests (max 20)
357      INTEGER  ::   ialert_id         ! number of the current alert
358      REAL(wp) ::   ztmelts           ! ice layer melting point
359      CHARACTER (len=30), DIMENSION(20) ::   cl_alname   ! name of alert
360      INTEGER           , DIMENSION(20) ::   inb_alp     ! number of alerts positive
361      !!-------------------------------------------------------------------
362
363      inb_altests = 10
364      inb_alp(:)  =  0
365
366      ! Alert if incompatible volume and concentration
367      ialert_id = 2 ! reference number of this alert
368      cl_alname(ialert_id) = ' Incompat vol and con         '    ! name of the alert
369      DO jl = 1, jpl
370         DO_2D_11_11
371            IF(  v_i(ji,jj,jl) /= 0._wp   .AND.   a_i(ji,jj,jl) == 0._wp   ) THEN
372               WRITE(numout,*) ' ALERTE 2 :   Incompatible volume and concentration '
373               inb_alp(ialert_id) = inb_alp(ialert_id) + 1
374            ENDIF
375         END_2D
376      END DO
377
378      ! Alerte if very thick ice
379      ialert_id = 3 ! reference number of this alert
380      cl_alname(ialert_id) = ' Very thick ice               ' ! name of the alert
381      jl = jpl 
382      DO_2D_11_11
383         IF(   h_i(ji,jj,jl)  >  50._wp   ) THEN
384            WRITE(numout,*) ' ALERTE 3 :   Very thick ice'
385            !CALL ice_prt( kt, ji, jj, 2, ' ALERTE 3 :   Very thick ice ' )
386            inb_alp(ialert_id) = inb_alp(ialert_id) + 1
387         ENDIF
388      END_2D
389
390      ! Alert if very fast ice
391      ialert_id = 4 ! reference number of this alert
392      cl_alname(ialert_id) = ' Very fast ice               ' ! name of the alert
393      DO_2D_11_11
394         IF(   MAX( ABS( u_ice(ji,jj) ), ABS( v_ice(ji,jj) ) ) > 2.  .AND.  &
395            &  at_i(ji,jj) > 0._wp   ) THEN
396            WRITE(numout,*) ' ALERTE 4 :   Very fast ice'
397            !CALL ice_prt( kt, ji, jj, 1, ' ALERTE 4 :   Very fast ice ' )
398            inb_alp(ialert_id) = inb_alp(ialert_id) + 1
399         ENDIF
400      END_2D
401
402      ! Alert on salt flux
403      ialert_id = 5 ! reference number of this alert
404      cl_alname(ialert_id) = ' High salt flux               ' ! name of the alert
405      DO_2D_11_11
406         IF( ABS( sfx (ji,jj) ) > 1.0e-2 ) THEN  ! = 1 psu/day for 1m ocean depth
407            WRITE(numout,*) ' ALERTE 5 :   High salt flux'
408            !CALL ice_prt( kt, ji, jj, 3, ' ALERTE 5 :   High salt flux ' )
409            inb_alp(ialert_id) = inb_alp(ialert_id) + 1
410         ENDIF
411      END_2D
412
413      ! Alert if there is ice on continents
414      ialert_id = 6 ! reference number of this alert
415      cl_alname(ialert_id) = ' Ice on continents           ' ! name of the alert
416      DO_2D_11_11
417         IF(   tmask(ji,jj,1) <= 0._wp   .AND.   at_i(ji,jj) > 0._wp   ) THEN
418            WRITE(numout,*) ' ALERTE 6 :   Ice on continents'
419            !CALL ice_prt( kt, ji, jj, 1, ' ALERTE 6 :   Ice on continents ' )
420            inb_alp(ialert_id) = inb_alp(ialert_id) + 1
421         ENDIF
422      END_2D
423
424!
425!     ! Alert if very fresh ice
426      ialert_id = 7 ! reference number of this alert
427      cl_alname(ialert_id) = ' Very fresh ice               ' ! name of the alert
428      DO jl = 1, jpl
429         DO_2D_11_11
430            IF( s_i(ji,jj,jl) < 0.1 .AND. a_i(ji,jj,jl) > 0._wp ) THEN
431               WRITE(numout,*) ' ALERTE 7 :   Very fresh ice'
432!                 CALL ice_prt(kt,ji,jj,1, ' ALERTE 7 :   Very fresh ice ' )
433               inb_alp(ialert_id) = inb_alp(ialert_id) + 1
434            ENDIF
435         END_2D
436      END DO
437!
438      ! Alert if qns very big
439      ialert_id = 8 ! reference number of this alert
440      cl_alname(ialert_id) = ' fnsolar very big             ' ! name of the alert
441      DO_2D_11_11
442         IF( ABS( qns(ji,jj) ) > 1500._wp .AND. at_i(ji,jj) > 0._wp ) THEN
443            !
444            WRITE(numout,*) ' ALERTE 8 :   Very high non-solar heat flux'
445            !CALL ice_prt( kt, ji, jj, 2, '   ')
446            inb_alp(ialert_id) = inb_alp(ialert_id) + 1
447            !
448         ENDIF
449      END_2D
450      !+++++
451
452!     ! Alert if too old ice
453      ialert_id = 9 ! reference number of this alert
454      cl_alname(ialert_id) = ' Very old   ice               ' ! name of the alert
455      DO jl = 1, jpl
456         DO_2D_11_11
457            IF ( ( ( ABS( o_i(ji,jj,jl) ) > rdt_ice ) .OR. &
458                   ( ABS( o_i(ji,jj,jl) ) < 0._wp) ) .AND. &
459                          ( a_i(ji,jj,jl) > 0._wp ) ) THEN
460               WRITE(numout,*) ' ALERTE 9 :   Wrong ice age'
461               !CALL ice_prt( kt, ji, jj, 1, ' ALERTE 9 :   Wrong ice age ')
462               inb_alp(ialert_id) = inb_alp(ialert_id) + 1
463            ENDIF
464         END_2D
465      END DO
466 
467      ! Alert if very warm ice
468      ialert_id = 10 ! reference number of this alert
469      cl_alname(ialert_id) = ' Very warm ice                ' ! name of the alert
470      inb_alp(ialert_id) = 0
471      DO jl = 1, jpl
472         DO_3D_11_11( 1, nlay_i )
473            ztmelts    =  -rTmlt * sz_i(ji,jj,jk,jl) + rt0
474            IF( t_i(ji,jj,jk,jl) > ztmelts  .AND.  v_i(ji,jj,jl) > 1.e-10   &
475               &                            .AND.  a_i(ji,jj,jl) > 0._wp   ) THEN
476               WRITE(numout,*) ' ALERTE 10 :   Very warm ice'
477              inb_alp(ialert_id) = inb_alp(ialert_id) + 1
478            ENDIF
479         END_3D
480      END DO
481
482      ! sum of the alerts on all processors
483      IF( lk_mpp ) THEN
484         DO ialert_id = 1, inb_altests
485            CALL mpp_sum('icectl', inb_alp(ialert_id))
486         END DO
487      ENDIF
488
489      ! print alerts
490      IF( lwp ) THEN
491         ialert_id = 1                                 ! reference number of this alert
492         cl_alname(ialert_id) = ' NO alerte 1      '   ! name of the alert
493         WRITE(numout,*) ' time step ',kt
494         WRITE(numout,*) ' All alerts at the end of ice model '
495         DO ialert_id = 1, inb_altests
496            WRITE(numout,*) ialert_id, cl_alname(ialert_id)//' : ', inb_alp(ialert_id), ' times ! '
497         END DO
498      ENDIF
499     !
500   END SUBROUTINE ice_ctl
501 
502   SUBROUTINE ice_prt( kt, ki, kj, kn, cd1 )
503      !!-------------------------------------------------------------------
504      !!                   ***  ROUTINE ice_prt ***
505      !!                 
506      !! ** Purpose :   Writes global ice state on the (i,j) point
507      !!                in ocean.ouput
508      !!                3 possibilities exist
509      !!                n = 1/-1 -> simple ice state
510      !!                n = 2    -> exhaustive state
511      !!                n = 3    -> ice/ocean salt fluxes
512      !!
513      !! ** input   :   point coordinates (i,j)
514      !!                n : number of the option
515      !!-------------------------------------------------------------------
516      INTEGER         , INTENT(in) ::   kt            ! ocean time step
517      INTEGER         , INTENT(in) ::   ki, kj, kn    ! ocean gridpoint indices
518      CHARACTER(len=*), INTENT(in) ::   cd1           !
519      !!
520      INTEGER :: jl, ji, jj
521      !!-------------------------------------------------------------------
522
523      DO ji = mi0(ki), mi1(ki)
524         DO jj = mj0(kj), mj1(kj)
525
526            WRITE(numout,*) ' time step ',kt,' ',cd1             ! print title
527
528            !----------------
529            !  Simple state
530            !----------------
531           
532            IF ( kn == 1 .OR. kn == -1 ) THEN
533               WRITE(numout,*) ' ice_prt - Point : ',ji,jj
534               WRITE(numout,*) ' ~~~~~~~~~~~~~~ '
535               WRITE(numout,*) ' Simple state '
536               WRITE(numout,*) ' masks s,u,v   : ', tmask(ji,jj,1), umask(ji,jj,1), vmask(ji,jj,1)
537               WRITE(numout,*) ' lat - long    : ', gphit(ji,jj), glamt(ji,jj)
538               WRITE(numout,*) ' - Ice drift   '
539               WRITE(numout,*) '   ~~~~~~~~~~~ '
540               WRITE(numout,*) ' u_ice(i-1,j)  : ', u_ice(ji-1,jj)
541               WRITE(numout,*) ' u_ice(i  ,j)  : ', u_ice(ji,jj)
542               WRITE(numout,*) ' v_ice(i  ,j-1): ', v_ice(ji,jj-1)
543               WRITE(numout,*) ' v_ice(i  ,j)  : ', v_ice(ji,jj)
544               WRITE(numout,*) ' strength      : ', strength(ji,jj)
545               WRITE(numout,*)
546               WRITE(numout,*) ' - Cell values '
547               WRITE(numout,*) '   ~~~~~~~~~~~ '
548               WRITE(numout,*) ' at_i          : ', at_i(ji,jj)       
549               WRITE(numout,*) ' ato_i         : ', ato_i(ji,jj)       
550               WRITE(numout,*) ' vt_i          : ', vt_i(ji,jj)       
551               WRITE(numout,*) ' vt_s          : ', vt_s(ji,jj)       
552               DO jl = 1, jpl
553                  WRITE(numout,*) ' - Category (', jl,')'
554                  WRITE(numout,*) ' a_i           : ', a_i(ji,jj,jl)
555                  WRITE(numout,*) ' h_i           : ', h_i(ji,jj,jl)
556                  WRITE(numout,*) ' h_s           : ', h_s(ji,jj,jl)
557                  WRITE(numout,*) ' v_i           : ', v_i(ji,jj,jl)
558                  WRITE(numout,*) ' v_s           : ', v_s(ji,jj,jl)
559                  WRITE(numout,*) ' e_s           : ', e_s(ji,jj,1:nlay_s,jl)
560                  WRITE(numout,*) ' e_i           : ', e_i(ji,jj,1:nlay_i,jl)
561                  WRITE(numout,*) ' t_su          : ', t_su(ji,jj,jl)
562                  WRITE(numout,*) ' t_snow        : ', t_s(ji,jj,1:nlay_s,jl)
563                  WRITE(numout,*) ' t_i           : ', t_i(ji,jj,1:nlay_i,jl)
564                  WRITE(numout,*) ' s_i           : ', s_i(ji,jj,jl)
565                  WRITE(numout,*) ' sv_i          : ', sv_i(ji,jj,jl)
566                  WRITE(numout,*)
567               END DO
568            ENDIF
569
570            !--------------------
571            !  Exhaustive state
572            !--------------------
573           
574            IF ( kn .EQ. 2 ) THEN
575               WRITE(numout,*) ' ice_prt - Point : ',ji,jj
576               WRITE(numout,*) ' ~~~~~~~~~~~~~~ '
577               WRITE(numout,*) ' Exhaustive state '
578               WRITE(numout,*) ' lat - long ', gphit(ji,jj), glamt(ji,jj)
579               WRITE(numout,*) 
580               WRITE(numout,*) ' - Cell values '
581               WRITE(numout,*) '   ~~~~~~~~~~~ '
582               WRITE(numout,*) ' at_i          : ', at_i(ji,jj)       
583               WRITE(numout,*) ' vt_i          : ', vt_i(ji,jj)       
584               WRITE(numout,*) ' vt_s          : ', vt_s(ji,jj)       
585               WRITE(numout,*) ' u_ice(i-1,j)  : ', u_ice(ji-1,jj)
586               WRITE(numout,*) ' u_ice(i  ,j)  : ', u_ice(ji,jj)
587               WRITE(numout,*) ' v_ice(i  ,j-1): ', v_ice(ji,jj-1)
588               WRITE(numout,*) ' v_ice(i  ,j)  : ', v_ice(ji,jj)
589               WRITE(numout,*) ' strength      : ', strength(ji,jj)
590               WRITE(numout,*) ' u_ice_b       : ', u_ice_b(ji,jj)    , ' v_ice_b       : ', v_ice_b(ji,jj) 
591               WRITE(numout,*)
592               
593               DO jl = 1, jpl
594                  WRITE(numout,*) ' - Category (',jl,')'
595                  WRITE(numout,*) '   ~~~~~~~~         ' 
596                  WRITE(numout,*) ' h_i        : ', h_i(ji,jj,jl)              , ' h_s        : ', h_s(ji,jj,jl)
597                  WRITE(numout,*) ' t_i        : ', t_i(ji,jj,1:nlay_i,jl)
598                  WRITE(numout,*) ' t_su       : ', t_su(ji,jj,jl)             , ' t_s        : ', t_s(ji,jj,1:nlay_s,jl)
599                  WRITE(numout,*) ' s_i        : ', s_i(ji,jj,jl)              , ' o_i        : ', o_i(ji,jj,jl)
600                  WRITE(numout,*) ' a_i        : ', a_i(ji,jj,jl)              , ' a_i_b      : ', a_i_b(ji,jj,jl)   
601                  WRITE(numout,*) ' v_i        : ', v_i(ji,jj,jl)              , ' v_i_b      : ', v_i_b(ji,jj,jl)   
602                  WRITE(numout,*) ' v_s        : ', v_s(ji,jj,jl)              , ' v_s_b      : ', v_s_b(ji,jj,jl) 
603                  WRITE(numout,*) ' e_i1       : ', e_i(ji,jj,1,jl)            , ' ei1        : ', e_i_b(ji,jj,1,jl) 
604                  WRITE(numout,*) ' e_i2       : ', e_i(ji,jj,2,jl)            , ' ei2_b      : ', e_i_b(ji,jj,2,jl) 
605                  WRITE(numout,*) ' e_snow     : ', e_s(ji,jj,1,jl)            , ' e_snow_b   : ', e_s_b(ji,jj,1,jl) 
606                  WRITE(numout,*) ' sv_i       : ', sv_i(ji,jj,jl)             , ' sv_i_b     : ', sv_i_b(ji,jj,jl)   
607                  WRITE(numout,*) ' oa_i       : ', oa_i(ji,jj,jl)             , ' oa_i_b     : ', oa_i_b(ji,jj,jl)
608               END DO !jl
609               
610               WRITE(numout,*)
611               WRITE(numout,*) ' - Heat / FW fluxes '
612               WRITE(numout,*) '   ~~~~~~~~~~~~~~~~ '
613               WRITE(numout,*) ' - Heat fluxes in and out the ice ***'
614               WRITE(numout,*) ' qsr_ini       : ', (1._wp-at_i_b(ji,jj)) * qsr(ji,jj) + SUM( a_i_b(ji,jj,:) * qsr_ice(ji,jj,:) )
615               WRITE(numout,*) ' qns_ini       : ', (1._wp-at_i_b(ji,jj)) * qns(ji,jj) + SUM( a_i_b(ji,jj,:) * qns_ice(ji,jj,:) )
616               WRITE(numout,*)
617               WRITE(numout,*) 
618               WRITE(numout,*) ' sst        : ', sst_m(ji,jj) 
619               WRITE(numout,*) ' sss        : ', sss_m(ji,jj) 
620               WRITE(numout,*) 
621               WRITE(numout,*) ' - Stresses '
622               WRITE(numout,*) '   ~~~~~~~~ '
623               WRITE(numout,*) ' utau_ice   : ', utau_ice(ji,jj) 
624               WRITE(numout,*) ' vtau_ice   : ', vtau_ice(ji,jj)
625               WRITE(numout,*) ' utau       : ', utau    (ji,jj) 
626               WRITE(numout,*) ' vtau       : ', vtau    (ji,jj)
627            ENDIF
628           
629            !---------------------
630            ! Salt / heat fluxes
631            !---------------------
632           
633            IF ( kn .EQ. 3 ) THEN
634               WRITE(numout,*) ' ice_prt - Point : ',ji,jj
635               WRITE(numout,*) ' ~~~~~~~~~~~~~~ '
636               WRITE(numout,*) ' - Salt / Heat Fluxes '
637               WRITE(numout,*) '   ~~~~~~~~~~~~~~~~ '
638               WRITE(numout,*) ' lat - long ', gphit(ji,jj), glamt(ji,jj)
639               WRITE(numout,*)
640               WRITE(numout,*) ' - Heat fluxes at bottom interface ***'
641               WRITE(numout,*) ' qsr       : ', qsr(ji,jj)
642               WRITE(numout,*) ' qns       : ', qns(ji,jj)
643               WRITE(numout,*)
644               WRITE(numout,*) ' hfx_mass     : ', hfx_thd(ji,jj) + hfx_dyn(ji,jj) + hfx_snw(ji,jj) + hfx_res(ji,jj)
645               WRITE(numout,*) ' qt_atm_oi    : ', qt_atm_oi(ji,jj)
646               WRITE(numout,*) ' qt_oce_ai    : ', qt_oce_ai(ji,jj)
647               WRITE(numout,*) ' dhc          : ', diag_heat(ji,jj)             
648               WRITE(numout,*)
649               WRITE(numout,*) ' hfx_dyn      : ', hfx_dyn(ji,jj)
650               WRITE(numout,*) ' hfx_thd      : ', hfx_thd(ji,jj)
651               WRITE(numout,*) ' hfx_res      : ', hfx_res(ji,jj)
652               WRITE(numout,*) ' qsb_ice_bot  : ', qsb_ice_bot(ji,jj) 
653               WRITE(numout,*) ' qlead        : ', qlead(ji,jj) * r1_rdtice
654               WRITE(numout,*)
655               WRITE(numout,*) ' - Salt fluxes at bottom interface ***'
656               WRITE(numout,*) ' emp       : ', emp    (ji,jj)
657               WRITE(numout,*) ' sfx       : ', sfx    (ji,jj)
658               WRITE(numout,*) ' sfx_res   : ', sfx_res(ji,jj)
659               WRITE(numout,*) ' sfx_bri   : ', sfx_bri(ji,jj)
660               WRITE(numout,*) ' sfx_dyn   : ', sfx_dyn(ji,jj)
661               WRITE(numout,*)
662               WRITE(numout,*) ' - Momentum fluxes '
663               WRITE(numout,*) ' utau      : ', utau(ji,jj) 
664               WRITE(numout,*) ' vtau      : ', vtau(ji,jj)
665            ENDIF
666            WRITE(numout,*) ' '
667            !
668         END DO
669      END DO
670      !
671   END SUBROUTINE ice_prt
672
673   SUBROUTINE ice_prt3D( cd_routine )
674      !!-------------------------------------------------------------------
675      !!                  ***  ROUTINE ice_prt3D ***
676      !!
677      !! ** Purpose : CTL prints of ice arrays in case sn_cfctl%prtctl is activated
678      !!
679      !!-------------------------------------------------------------------
680      CHARACTER(len=*), INTENT(in) ::   cd_routine    ! name of the routine
681      INTEGER                      ::   jk, jl        ! dummy loop indices
682     
683      CALL prt_ctl_info(' ========== ')
684      CALL prt_ctl_info( cd_routine )
685      CALL prt_ctl_info(' ========== ')
686      CALL prt_ctl_info(' - Cell values : ')
687      CALL prt_ctl_info('   ~~~~~~~~~~~~~ ')
688      CALL prt_ctl(tab2d_1=e1e2t      , clinfo1=' cell area   :')
689      CALL prt_ctl(tab2d_1=at_i       , clinfo1=' at_i        :')
690      CALL prt_ctl(tab2d_1=ato_i      , clinfo1=' ato_i       :')
691      CALL prt_ctl(tab2d_1=vt_i       , clinfo1=' vt_i        :')
692      CALL prt_ctl(tab2d_1=vt_s       , clinfo1=' vt_s        :')
693      CALL prt_ctl(tab2d_1=divu_i     , clinfo1=' divu_i      :')
694      CALL prt_ctl(tab2d_1=delta_i    , clinfo1=' delta_i     :')
695      CALL prt_ctl(tab2d_1=stress1_i  , clinfo1=' stress1_i   :')
696      CALL prt_ctl(tab2d_1=stress2_i  , clinfo1=' stress2_i   :')
697      CALL prt_ctl(tab2d_1=stress12_i , clinfo1=' stress12_i  :')
698      CALL prt_ctl(tab2d_1=strength   , clinfo1=' strength    :')
699      CALL prt_ctl(tab2d_1=delta_i    , clinfo1=' delta_i     :')
700      CALL prt_ctl(tab2d_1=u_ice      , clinfo1=' u_ice       :', tab2d_2=v_ice      , clinfo2=' v_ice       :')
701       
702      DO jl = 1, jpl
703         CALL prt_ctl_info(' ')
704         CALL prt_ctl_info(' - Category : ', ivar1=jl)
705         CALL prt_ctl_info('   ~~~~~~~~~~')
706         CALL prt_ctl(tab2d_1=h_i        (:,:,jl)        , clinfo1= ' h_i         : ')
707         CALL prt_ctl(tab2d_1=h_s        (:,:,jl)        , clinfo1= ' h_s         : ')
708         CALL prt_ctl(tab2d_1=t_su       (:,:,jl)        , clinfo1= ' t_su        : ')
709         CALL prt_ctl(tab2d_1=t_s        (:,:,1,jl)      , clinfo1= ' t_snow      : ')
710         CALL prt_ctl(tab2d_1=s_i        (:,:,jl)        , clinfo1= ' s_i         : ')
711         CALL prt_ctl(tab2d_1=o_i        (:,:,jl)        , clinfo1= ' o_i         : ')
712         CALL prt_ctl(tab2d_1=a_i        (:,:,jl)        , clinfo1= ' a_i         : ')
713         CALL prt_ctl(tab2d_1=v_i        (:,:,jl)        , clinfo1= ' v_i         : ')
714         CALL prt_ctl(tab2d_1=v_s        (:,:,jl)        , clinfo1= ' v_s         : ')
715         CALL prt_ctl(tab2d_1=e_i        (:,:,1,jl)      , clinfo1= ' e_i1        : ')
716         CALL prt_ctl(tab2d_1=e_s        (:,:,1,jl)      , clinfo1= ' e_snow      : ')
717         CALL prt_ctl(tab2d_1=sv_i       (:,:,jl)        , clinfo1= ' sv_i        : ')
718         CALL prt_ctl(tab2d_1=oa_i       (:,:,jl)        , clinfo1= ' oa_i        : ')
719         
720         DO jk = 1, nlay_i
721            CALL prt_ctl_info(' - Layer : ', ivar1=jk)
722            CALL prt_ctl(tab2d_1=t_i(:,:,jk,jl) , clinfo1= ' t_i       : ')
723         END DO
724      END DO
725     
726      CALL prt_ctl_info(' ')
727      CALL prt_ctl_info(' - Heat / FW fluxes : ')
728      CALL prt_ctl_info('   ~~~~~~~~~~~~~~~~~~ ')
729      CALL prt_ctl(tab2d_1=sst_m  , clinfo1= ' sst   : ', tab2d_2=sss_m     , clinfo2= ' sss       : ')
730      CALL prt_ctl(tab2d_1=qsr    , clinfo1= ' qsr   : ', tab2d_2=qns       , clinfo2= ' qns       : ')
731      CALL prt_ctl(tab2d_1=emp    , clinfo1= ' emp   : ', tab2d_2=sfx       , clinfo2= ' sfx       : ')
732     
733      CALL prt_ctl_info(' ')
734      CALL prt_ctl_info(' - Stresses : ')
735      CALL prt_ctl_info('   ~~~~~~~~~~ ')
736      CALL prt_ctl(tab2d_1=utau       , clinfo1= ' utau      : ', tab2d_2=vtau       , clinfo2= ' vtau      : ')
737      CALL prt_ctl(tab2d_1=utau_ice   , clinfo1= ' utau_ice  : ', tab2d_2=vtau_ice   , clinfo2= ' vtau_ice  : ')
738     
739   END SUBROUTINE ice_prt3D
740     
741#else
742   !!----------------------------------------------------------------------
743   !!   Default option         Empty Module           No SI3 sea-ice model
744   !!----------------------------------------------------------------------
745#endif
746
747   !!======================================================================
748END MODULE icectl
Note: See TracBrowser for help on using the repository browser.