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.
icestp.F90 in NEMO/branches/UKMO/dev_r9950_GO8_package/src/ICE – NEMO

source: NEMO/branches/UKMO/dev_r9950_GO8_package/src/ICE/icestp.F90 @ 10322

Last change on this file since 10322 was 10322, checked in by davestorkey, 5 years ago

UKMO/dev_r9950_GO8_package: Update to be relative to rev 10321 of NEMO4_beta_mirror branch.

  • Property svn:keywords set to Id
File size: 22.7 KB
Line 
1MODULE icestp
2   !!======================================================================
3   !!                       ***  MODULE  icestp  ***
4   !! sea ice : Master routine for all the sea ice model
5   !!=====================================================================
6   !!
7   !! The sea ice model SI3 (Sea Ice modelling Integrated Initiative),
8   !!                        aka Sea Ice cube for its nickname
9   !!
10   !!    is originally based on LIM3, developed in Louvain-la-Neuve by:
11   !!       * Martin Vancoppenolle (UCL-ASTR, Belgium)
12   !!       * Sylvain Bouillon (UCL-ASTR, Belgium)
13   !!       * Miguel Angel Morales Maqueda (NOC-L, UK)
14   !!      thanks to valuable earlier work by
15   !!       * Thierry Fichefet
16   !!       * Hugues Goosse
17   !!      thanks also to the following persons who contributed
18   !!       * Gurvan Madec, Claude Talandier, Christian Ethe (LOCEAN, France)
19   !!       * Xavier Fettweis (UCL-ASTR), Ralph Timmermann (AWI, Germany)
20   !!       * Bill Lipscomb (LANL), Cecilia Bitz (UWa) and Elisabeth Hunke (LANL), USA.
21   !!
22   !! SI3 has been made possible by a handful of persons who met as working group
23   !!      (from France, Belgium, UK and Italy)
24   !!    * Clement Rousset, Martin Vancoppenolle & Gurvan Madec (LOCEAN, France)
25   !!    * Matthieu Chevalier & David Salas (Meteo France, France)
26   !!    * Gilles Garric (Mercator Ocean, France)
27   !!    * Thierry Fichefet & Francois Massonnet (UCL, Belgium)
28   !!    * Ed Blockley & Jeff Ridley (Met Office, UK)
29   !!    * Danny Feltham & David Schroeder (CPOM, UK)
30   !!    * Yevgeny Aksenov (NOC, UK)
31   !!    * Paul Holland (BAS, UK)
32   !!    * Dorotea Iovino (CMCC, Italy)
33   !!======================================================================
34   !! History :  4.0  !  2018     (C. Rousset)      Original code SI3
35   !!----------------------------------------------------------------------
36#if defined key_si3
37   !!----------------------------------------------------------------------
38   !!   'key_si3'                                       SI3 sea-ice model
39   !!----------------------------------------------------------------------
40   !!   ice_stp       : sea-ice model time-stepping and update ocean SBC over ice-covered area
41   !!   ice_init      : initialize sea-ice
42   !!----------------------------------------------------------------------
43   USE oce            ! ocean dynamics and tracers
44   USE dom_oce        ! ocean space and time domain
45   USE c1d            ! 1D vertical configuration
46   USE ice            ! sea-ice: variables
47   USE ice1D          ! sea-ice: thermodynamical 1D variables
48   !
49   USE phycst         ! Define parameters for the routines
50   USE eosbn2         ! equation of state
51   USE sbc_oce        ! Surface boundary condition: ocean fields
52   USE sbc_ice        ! Surface boundary condition: ice   fields
53   !
54   USE iceforcing     ! sea-ice: Surface boundary condition       !!gm why not icesbc module name
55   USE icedyn         ! sea-ice: dynamics
56   USE icethd         ! sea-ice: thermodynamics
57   USE icecor         ! sea-ice: corrections
58   USE iceupdate      ! sea-ice: sea surface boundary condition update
59   USE icedia         ! sea-ice: budget diagnostics
60   USE icewri         ! sea-ice: outputs
61   USE icerst         ! sea-ice: restarts
62   USE icevar         ! sea-ice: operations
63   USE icectl         ! sea-ice: control
64   USE iceistate      ! sea-ice: initial state
65   USE iceitd         ! sea-ice: remapping thickness distribution
66   USE icealb         ! sea-ice: albedo
67   !
68   USE bdy_oce , ONLY : ln_bdy   ! flag for bdy
69   USE bdyice         ! unstructured open boundary data for sea-ice
70# if defined key_agrif
71   USE agrif_ice
72   USE agrif_ice_interp
73# endif
74   !
75   USE in_out_manager ! I/O manager
76   USE iom            ! I/O manager library
77   USE lib_mpp        ! MPP library
78   USE lib_fortran    ! fortran utilities (glob_sum + no signed zero)
79   USE timing         ! Timing
80   USE prtctl         ! Print control
81
82   IMPLICIT NONE
83   PRIVATE
84
85   PUBLIC   ice_stp    ! called by sbcmod.F90
86   PUBLIC   ice_init   ! called by sbcmod.F90
87
88   !! * Substitutions
89#  include "vectopt_loop_substitute.h90"
90   !!----------------------------------------------------------------------
91   !! NEMO/ICE 4.0 , NEMO Consortium (2018)
92   !! $Id$
93   !! Software governed by the CeCILL license (see ./LICENSE)
94   !!----------------------------------------------------------------------
95CONTAINS
96
97   SUBROUTINE ice_stp( kt, ksbc )
98      !!---------------------------------------------------------------------
99      !!                  ***  ROUTINE ice_stp  ***
100      !!
101      !! ** Purpose :   sea-ice model time-stepping and update ocean surface
102      !!              boundary condition over ice-covered area
103      !!
104      !! ** Method  :   ice model time stepping
105      !!              - call the ice dynamics routine
106      !!              - call the ice advection/diffusion routine
107      !!              - call the ice thermodynamics routine
108      !!              - call the routine that computes mass and
109      !!                heat fluxes at the ice/ocean interface
110      !!              - save the outputs
111      !!              - save the outputs for restart when necessary
112      !!
113      !! ** Action  : - time evolution of the LIM sea-ice model
114      !!              - update all sbc variables below sea-ice:
115      !!                utau, vtau, taum, wndm, qns , qsr, emp , sfx
116      !!---------------------------------------------------------------------
117      INTEGER, INTENT(in) ::   kt      ! ocean time step
118      INTEGER, INTENT(in) ::   ksbc    ! flux formulation (user defined, bulk, or Pure Coupled)
119      !
120      INTEGER ::   jl   ! dummy loop index
121      !!----------------------------------------------------------------------
122      !
123      IF( ln_timing )   CALL timing_start('ice_stp')
124      !
125      !                                      !-----------------------!
126      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN   ! --- Ice time step --- !
127         !                                   !-----------------------!
128         !
129         kt_ice = kt                              ! -- Ice model time step
130         !
131         u_oce(:,:) = ssu_m(:,:)                  ! -- mean surface ocean current
132         v_oce(:,:) = ssv_m(:,:)
133         !
134         CALL eos_fzp( sss_m(:,:) , t_bo(:,:) )   ! -- freezing temperature [Kelvin] (set to rt0 over land)
135         t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1) + rt0 * ( 1._wp - tmask(:,:,1) )
136         !
137         !                          !==  AGRIF Parent to Child  ==!
138#if defined key_agrif
139         !                              ! nbstep_ice ranges from 1 to the nb of child ocean steps inside one parent ice step
140         IF( .NOT. Agrif_Root() )       nbstep_ice = MOD( nbstep_ice, Agrif_irhot() * Agrif_Parent(nn_fsbc) / nn_fsbc ) + 1
141         !                              ! these calls must remain here for restartability purposes
142                                        CALL agrif_interp_ice( 'T' ) 
143                                        CALL agrif_interp_ice( 'U' )
144                                        CALL agrif_interp_ice( 'V' )
145#endif
146                                        CALL store_fields             ! Store now ice values
147         !
148         !------------------------------------------------!
149         ! --- Dynamical coupling with the atmosphere --- !
150         !------------------------------------------------!
151         ! It provides the following fields used in sea ice model:
152         !    utau_ice, vtau_ice = surface ice stress [N/m2]
153         !------------------------------------------------!
154                                        CALL ice_forcing_tau( kt, ksbc, utau_ice, vtau_ice )         
155         !-------------------------------------!
156         ! --- ice dynamics and advection  --- !
157         !-------------------------------------!
158                                        CALL diag_set0                ! set diag of mass, heat and salt fluxes to 0
159                                        CALL ice_rst_opn( kt )        ! Open Ice restart file (if necessary)
160         !
161         IF( ln_icedyn .AND. .NOT.lk_c1d )   &
162            &                           CALL ice_dyn( kt )            ! -- Ice dynamics
163         !
164         !                          !==  lateral boundary conditions  ==!
165         IF( ln_icethd .AND. ln_bdy )   CALL bdy_ice( kt )            ! -- bdy ice thermo
166         !
167         !                          !==  previous lead fraction and ice volume for flux calculations
168                                        CALL ice_var_glo2eqv          ! h_i and h_s for ice albedo calculation
169                                        CALL ice_var_agg(1)           ! at_i for coupling
170                                        CALL store_fields             ! Store now ice values
171         !
172         !------------------------------------------------------!
173         ! --- Thermodynamical coupling with the atmosphere --- !
174         !------------------------------------------------------!
175         ! It provides the following fields used in sea ice model:
176         !    emp_oce , emp_ice    = E-P over ocean and sea ice                    [Kg/m2/s]
177         !    sprecip              = solid precipitation                           [Kg/m2/s]
178         !    evap_ice             = sublimation                                   [Kg/m2/s]
179         !    qsr_tot , qns_tot    = solar & non solar heat flux (total)           [W/m2]
180         !    qsr_ice , qns_ice    = solar & non solar heat flux over ice          [W/m2]
181         !    dqns_ice             = non solar  heat sensistivity                  [W/m2]
182         !    qemp_oce, qemp_ice,  = sensible heat (associated with evap & precip) [W/m2]
183         !    qprec_ice, qevap_ice
184         !------------------------------------------------------!
185                                        CALL ice_forcing_flx( kt, ksbc )
186         !----------------------------!
187         ! --- ice thermodynamics --- !
188         !----------------------------!
189         IF( ln_icethd )                CALL ice_thd( kt )            ! -- Ice thermodynamics     
190         !
191         IF( ln_icethd )                CALL ice_cor( kt , 2 )        ! -- Corrections
192         !
193                                        CALL ice_var_glo2eqv          ! necessary calls (at least for coupling)
194                                        CALL ice_var_agg( 2 )         ! necessary calls (at least for coupling)
195         !
196!! clem: one should switch the calculation of the fluxes onto the parent grid but the following calls do not work
197!!       moreover it should only be called at the update frequency (as in agrif_ice_update.F90)
198!# if defined key_agrif
199!         IF( .NOT. Agrif_Root() )     CALL Agrif_ChildGrid_To_ParentGrid()
200!# endif
201                                        CALL ice_update_flx( kt )     ! -- Update ocean surface mass, heat and salt fluxes
202!# if defined key_agrif
203!         IF( .NOT. Agrif_Root() )     CALL Agrif_ParentGrid_To_ChildGrid()
204!# endif
205         IF( ln_icediahsb )             CALL ice_dia( kt )            ! -- Diagnostics outputs
206         !
207                                        CALL ice_wri( kt )            ! -- Ice outputs
208         !
209         IF( lrst_ice )                 CALL ice_rst_write( kt )      ! -- Ice restart file
210         !
211         IF( ln_icectl )                CALL ice_ctl( kt )            ! -- alerts in case of model crash
212         !
213      ENDIF   ! End sea-ice time step only
214
215      !-------------------------!
216      ! --- Ocean time step --- !
217      !-------------------------!
218      IF( ln_icedyn )                   CALL ice_update_tau( kt, ub(:,:,1), vb(:,:,1) )   ! -- update surface ocean stresses
219!!gm   remark, the ocean-ice stress is not saved in ice diag call above .....  find a solution!!!
220      !
221      IF( ln_timing )   CALL timing_stop('ice_stp')
222      !
223   END SUBROUTINE ice_stp
224
225
226   SUBROUTINE ice_init
227      !!----------------------------------------------------------------------
228      !!                  ***  ROUTINE ice_init  ***
229      !!
230      !! ** purpose :   Initialize sea-ice parameters
231      !!----------------------------------------------------------------------
232      INTEGER :: ji, jj, ierr
233      !!----------------------------------------------------------------------
234      IF(lwp) WRITE(numout,*)
235      IF(lwp) WRITE(numout,*) 'ice_init: Arrays allocation & Initialization off all routines & init state' 
236      IF(lwp) WRITE(numout,*) '~~~~~~~~'
237      !
238      !                                ! Open the reference and configuration namelist files and namelist output file
239      CALL ctl_opn( numnam_ice_ref, 'namelist_ice_ref',    'OLD',     'FORMATTED', 'SEQUENTIAL', -1, numout, lwp )
240      CALL ctl_opn( numnam_ice_cfg, 'namelist_ice_cfg',    'OLD',     'FORMATTED', 'SEQUENTIAL', -1, numout, lwp )
241      IF(lwm) CALL ctl_opn( numoni, 'output.namelist.ice', 'UNKNOWN', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp, 1 )
242      !
243      CALL par_init                ! set some ice run parameters
244      !
245      !                                ! Allocate the ice arrays (sbc_ice already allocated in sbc_init)
246      ierr =        ice_alloc        ()      ! ice variables
247      ierr = ierr + sbc_ice_alloc    ()      ! surface forcing
248      ierr = ierr + ice1D_alloc      ()      ! thermodynamics
249      !
250      IF( lk_mpp    )   CALL mpp_sum( ierr )
251      IF( ierr /= 0 )   CALL ctl_stop('STOP', 'ice_init : unable to allocate ice arrays')
252      !
253      CALL ice_itd_init                ! ice thickness distribution initialization
254      !
255      CALL ice_thd_init                ! set ice thermodynics parameters (clem: important to call it first for melt ponds)
256      !
257      !                                ! Initial sea-ice state
258      IF( .NOT. ln_rstart ) THEN              ! start from rest: sea-ice deduced from sst
259         CALL ice_istate_init
260         CALL ice_istate
261      ELSE                                    ! start from a restart file
262         CALL ice_rst_read
263      ENDIF
264      CALL ice_var_glo2eqv
265      CALL ice_var_agg(1)
266      !
267      CALL ice_forcing_init            ! set ice-ocean and ice-atm. coupling parameters
268      !
269      CALL ice_dyn_init                ! set ice dynamics parameters
270      !
271      CALL ice_update_init             ! ice surface boundary condition
272      !
273      CALL ice_alb_init                ! ice surface albedo
274      !
275      CALL ice_dia_init                ! initialization for diags
276      !
277      fr_i  (:,:)   = at_i(:,:)        ! initialisation of sea-ice fraction
278      tn_ice(:,:,:) = t_su(:,:,:)      ! initialisation of surface temp for coupled simu
279      !
280      !                                ! set max concentration in both hemispheres
281      WHERE( gphit(:,:) > 0._wp )   ;   rn_amax_2d(:,:) = rn_amax_n  ! NH
282      ELSEWHERE                     ;   rn_amax_2d(:,:) = rn_amax_s  ! SH
283      END WHERE
284
285      IF( ln_rstart )   CALL iom_close( numrir )  ! close input ice restart file
286      !
287   END SUBROUTINE ice_init
288
289
290   SUBROUTINE par_init
291      !!-------------------------------------------------------------------
292      !!                  ***  ROUTINE par_init ***
293      !!
294      !! ** Purpose :   Definition generic parameters for ice model
295      !!
296      !! ** Method  :   Read namelist and check the parameter
297      !!                values called at the first timestep (nit000)
298      !!
299      !! ** input   :   Namelist nampar
300      !!-------------------------------------------------------------------
301      INTEGER  ::   ios                 ! Local integer
302      !!
303      NAMELIST/nampar/ jpl, nlay_i, nlay_s, nn_virtual_itd, ln_icedyn, ln_icethd, rn_amax_n, rn_amax_s,  &
304         &             cn_icerst_in, cn_icerst_indir, cn_icerst_out, cn_icerst_outdir
305      !!-------------------------------------------------------------------
306      !
307      REWIND( numnam_ice_ref )      ! Namelist nampar in reference namelist : Parameters for ice
308      READ  ( numnam_ice_ref, nampar, IOSTAT = ios, ERR = 901)
309901   IF( ios /= 0 )   CALL ctl_nam ( ios , 'nampar in reference namelist', lwp )
310      REWIND( numnam_ice_cfg )      ! Namelist nampar in configuration namelist : Parameters for ice
311      READ  ( numnam_ice_cfg, nampar, IOSTAT = ios, ERR = 902 )
312902   IF( ios > 0 )   CALL ctl_nam ( ios , 'nampar in configuration namelist', lwp )
313      IF(lwm) WRITE( numoni, nampar )
314      !
315      IF(lwp) THEN                  ! control print
316         WRITE(numout,*)
317         WRITE(numout,*) '   par_init: ice parameters shared among all the routines'
318         WRITE(numout,*) '   ~~~~~~~~'
319         WRITE(numout,*) '      Namelist nampar: '
320         WRITE(numout,*) '         number of ice  categories                           jpl       = ', jpl
321         WRITE(numout,*) '         number of ice  layers                               nlay_i    = ', nlay_i
322         WRITE(numout,*) '         number of snow layers                               nlay_s    = ', nlay_s
323         WRITE(numout,*) '         virtual ITD param for jpl=1 (1-3) or not (0)   nn_virtual_itd = ', nn_virtual_itd
324         WRITE(numout,*) '         Ice dynamics       (T) or not (F)                   ln_icedyn = ', ln_icedyn
325         WRITE(numout,*) '         Ice thermodynamics (T) or not (F)                   ln_icethd = ', ln_icethd
326         WRITE(numout,*) '         maximum ice concentration for NH                              = ', rn_amax_n 
327         WRITE(numout,*) '         maximum ice concentration for SH                              = ', rn_amax_s
328      ENDIF
329      !                                        !--- check consistency
330      IF ( jpl > 1 .AND. nn_virtual_itd == 1 ) THEN
331         nn_virtual_itd = 0
332         IF(lwp) WRITE(numout,*)
333         IF(lwp) WRITE(numout,*) '   nn_virtual_itd forced to 0 as jpl>1, no need with multiple categories to emulate them'
334      ENDIF
335      !
336      IF( ln_cpl .AND. nn_cats_cpl /= 1 .AND. nn_cats_cpl /= jpl ) THEN
337         CALL ctl_stop( 'STOP', 'par_init: in coupled mode, nn_cats_cpl should be either 1 or jpl' )
338      ENDIF
339      !
340      IF( ln_bdy .AND. ln_icediachk )   CALL ctl_warn('par_init: online conservation check does not work with BDY')
341      !
342      rdt_ice   = REAL(nn_fsbc) * rdt          !--- sea-ice timestep and its inverse
343      r1_rdtice = 1._wp / rdt_ice
344      IF(lwp) WRITE(numout,*)
345      IF(lwp) WRITE(numout,*) '      ice timestep rdt_ice = nn_fsbc*rdt = ', rdt_ice
346      !
347      r1_nlay_i = 1._wp / REAL( nlay_i, wp )   !--- inverse of nlay_i and nlay_s
348      r1_nlay_s = 1._wp / REAL( nlay_s, wp )
349      !
350   END SUBROUTINE par_init
351
352
353   SUBROUTINE store_fields
354      !!----------------------------------------------------------------------
355      !!                  ***  ROUTINE store_fields  ***
356      !!
357      !! ** purpose :  store ice variables at "before" time step
358      !!----------------------------------------------------------------------
359      INTEGER  ::   ji, jj, jl      ! dummy loop index
360      !!----------------------------------------------------------------------
361      !
362      a_i_b (:,:,:)   = a_i (:,:,:)     ! ice area
363      v_i_b (:,:,:)   = v_i (:,:,:)     ! ice volume
364      v_s_b (:,:,:)   = v_s (:,:,:)     ! snow volume
365      sv_i_b(:,:,:)   = sv_i(:,:,:)     ! salt content
366      oa_i_b(:,:,:)   = oa_i(:,:,:)     ! areal age content
367      e_s_b (:,:,:,:) = e_s (:,:,:,:)   ! snow thermal energy
368      e_i_b (:,:,:,:) = e_i (:,:,:,:)   ! ice thermal energy
369      WHERE( a_i_b(:,:,:) >= epsi20 )
370         h_i_b(:,:,:) = v_i_b (:,:,:) / a_i_b(:,:,:)   ! ice thickness
371         h_s_b(:,:,:) = v_s_b (:,:,:) / a_i_b(:,:,:)   ! snw thickness
372      ELSEWHERE
373         h_i_b(:,:,:) = 0._wp
374         h_s_b(:,:,:) = 0._wp
375      END WHERE
376      !
377      ! ice velocities & total concentration
378      at_i_b(:,:)  = SUM( a_i_b(:,:,:), dim=3 )
379      u_ice_b(:,:) = u_ice(:,:)
380      v_ice_b(:,:) = v_ice(:,:)
381      !
382   END SUBROUTINE store_fields
383
384
385   SUBROUTINE diag_set0
386      !!----------------------------------------------------------------------
387      !!                  ***  ROUTINE diag_set0  ***
388      !!
389      !! ** purpose :  set ice-ocean and ice-atm. fluxes to zeros at the beggining
390      !!               of the time step
391      !!----------------------------------------------------------------------
392      INTEGER  ::   ji, jj      ! dummy loop index
393      !!----------------------------------------------------------------------
394      sfx    (:,:) = 0._wp   ;
395      sfx_bri(:,:) = 0._wp   ;   sfx_lam(:,:) = 0._wp
396      sfx_sni(:,:) = 0._wp   ;   sfx_opw(:,:) = 0._wp
397      sfx_bog(:,:) = 0._wp   ;   sfx_dyn(:,:) = 0._wp
398      sfx_bom(:,:) = 0._wp   ;   sfx_sum(:,:) = 0._wp
399      sfx_res(:,:) = 0._wp   ;   sfx_sub(:,:) = 0._wp
400      !
401      wfx_snw(:,:) = 0._wp   ;   wfx_ice(:,:) = 0._wp
402      wfx_sni(:,:) = 0._wp   ;   wfx_opw(:,:) = 0._wp
403      wfx_bog(:,:) = 0._wp   ;   wfx_dyn(:,:) = 0._wp
404      wfx_bom(:,:) = 0._wp   ;   wfx_sum(:,:) = 0._wp
405      wfx_res(:,:) = 0._wp   ;   wfx_sub(:,:) = 0._wp
406      wfx_spr(:,:) = 0._wp   ;   wfx_lam(:,:) = 0._wp 
407      wfx_snw_dyn(:,:) = 0._wp ; wfx_snw_sum(:,:) = 0._wp
408      wfx_snw_sub(:,:) = 0._wp ; wfx_ice_sub(:,:) = 0._wp
409      wfx_snw_sni(:,:) = 0._wp 
410      wfx_pnd(:,:) = 0._wp
411
412      hfx_thd(:,:) = 0._wp   ;
413      hfx_snw(:,:) = 0._wp   ;   hfx_opw(:,:) = 0._wp
414      hfx_bog(:,:) = 0._wp   ;   hfx_dyn(:,:) = 0._wp
415      hfx_bom(:,:) = 0._wp   ;   hfx_sum(:,:) = 0._wp
416      hfx_res(:,:) = 0._wp   ;   hfx_sub(:,:) = 0._wp
417      hfx_spr(:,:) = 0._wp   ;   hfx_dif(:,:) = 0._wp
418      hfx_err_rem(:,:) = 0._wp
419      hfx_err_dif(:,:) = 0._wp
420      wfx_err_sub(:,:) = 0._wp
421      !
422      afx_tot(:,:) = 0._wp   ;
423      !
424      diag_heat(:,:) = 0._wp ;   diag_sice(:,:) = 0._wp
425      diag_vice(:,:) = 0._wp ;   diag_vsnw(:,:) = 0._wp
426
427      ! SIMIP diagnostics
428      qcn_ice_bot(:,:,:) = 0._wp ; qcn_ice_top(:,:,:) = 0._wp ! conductive fluxes
429      t_si       (:,:,:) = rt0   ! temp at the ice-snow interface
430
431      tau_icebfr(:,:)   = 0._wp   ! landfast ice param only (clem: important to keep the init here)
432      cnd_ice   (:,:,:) = 0._wp   ! initialisation: effective conductivity at the top of ice/snow (Jules coupling)
433      qtr_ice_bot(:,:,:) = 0._wp  ! initialization: part of solar radiation transmitted through the ice needed at least for outputs
434      !
435      ! for control checks (ln_icediachk)
436      diag_trp_vi(:,:) = 0._wp   ;   diag_trp_vs(:,:) = 0._wp
437      diag_trp_ei(:,:) = 0._wp   ;   diag_trp_es(:,:) = 0._wp
438      diag_trp_sv(:,:) = 0._wp
439
440   END SUBROUTINE diag_set0
441
442#else
443   !!----------------------------------------------------------------------
444   !!   Default option           Dummy module         NO SI3 sea-ice model
445   !!----------------------------------------------------------------------
446CONTAINS
447   SUBROUTINE ice_stp ( kt, ksbc )     ! Dummy routine
448      INTEGER, INTENT(in) ::   kt, ksbc
449      WRITE(*,*) 'ice_stp: You should not have seen this print! error?', kt
450   END SUBROUTINE ice_stp
451   SUBROUTINE ice_init                 ! Dummy routine
452   END SUBROUTINE ice_init
453#endif
454
455   !!======================================================================
456END MODULE icestp
Note: See TracBrowser for help on using the repository browser.