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.
limdyn.F90 in branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2016/dev_v3_6_STABLE_r6506_AGRIF_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90 @ 6746

Last change on this file since 6746 was 6746, checked in by clem, 8 years ago

landfast ice parameterization + update from trunk + removing useless dom_ice.F90 and limmsh.F90 and limwri_dimg.h90
  • Property svn:keywords set to Id
File size: 9.8 KB
Line 
1MODULE limdyn
2   !!======================================================================
3   !!                     ***  MODULE  limdyn  ***
4   !!   Sea-Ice dynamics : 
5   !!======================================================================
6   !! history :  1.0  ! 2002-08  (C. Ethe, G. Madec)  original VP code
7   !!            3.0  ! 2007-03  (MA Morales Maqueda, S. Bouillon, M. Vancoppenolle)  LIM3: EVP-Cgrid
8   !!            3.5  ! 2011-02  (G. Madec) dynamical allocation
9   !!----------------------------------------------------------------------
10#if defined key_lim3
11   !!----------------------------------------------------------------------
12   !!   'key_lim3' :                                 LIM3 sea-ice model
13   !!----------------------------------------------------------------------
14   !!    lim_dyn      : computes ice velocities
15   !!    lim_dyn_init : initialization and namelist read
16   !!----------------------------------------------------------------------
17   USE phycst           ! physical constants
18   USE dom_oce          ! ocean space and time domain
19   USE sbc_ice          ! Surface boundary condition: ice   fields
20   USE ice              ! LIM-3 variables
21   USE limrhg           ! LIM-3 rheology
22   USE lbclnk           ! lateral boundary conditions - MPP exchanges
23   USE lib_mpp          ! MPP library
24   USE wrk_nemo         ! work arrays
25   USE in_out_manager   ! I/O manager
26   USE prtctl           ! Print control
27   USE lib_fortran      ! glob_sum
28   USE timing          ! Timing
29   USE limcons        ! conservation tests
30   USE limvar
31
32   IMPLICIT NONE
33   PRIVATE
34
35   PUBLIC   lim_dyn   ! routine called by ice_step
36
37   !! * Substitutions
38#  include "vectopt_loop_substitute.h90"
39   !!----------------------------------------------------------------------
40   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
41   !! $Id$
42   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
43   !!----------------------------------------------------------------------
44CONTAINS
45
46   SUBROUTINE lim_dyn( kt )
47      !!-------------------------------------------------------------------
48      !!               ***  ROUTINE lim_dyn  ***
49      !!               
50      !! ** Purpose :   compute ice velocity
51      !!               
52      !! ** Method  :
53      !!
54      !! ** Action  : - Initialisation
55      !!              - Call of the dynamic routine for each hemisphere
56      !!------------------------------------------------------------------------------------
57      INTEGER, INTENT(in) ::   kt     ! number of iteration
58      !!
59      INTEGER  :: jl, jk ! dummy loop indices
60      REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 
61     !!---------------------------------------------------------------------
62
63      IF( nn_timing == 1 )  CALL timing_start('limdyn')
64
65      CALL lim_var_agg(1)                      ! aggregate ice categories
66
67      IF( kt == nit000 )   CALL lim_dyn_init   ! Initialization (first time-step only)
68      !
69      ! conservation test
70      IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limdyn', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b)
71     
72      ! ice velocities before rheology
73      u_ice_b(:,:) = u_ice(:,:) * umask(:,:,1)
74      v_ice_b(:,:) = v_ice(:,:) * vmask(:,:,1)
75     
76      ! Landfast ice parameterization: define max bottom friction
77      tau_icebfr(:,:) = 0._wp
78      IF( ln_landfast ) THEN
79         DO jl = 1, jpl
80            WHERE( ht_i(:,:,jl) > bathy(:,:) * rn_gamma )  tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr
81         END DO
82      ENDIF
83     
84      ! Rheology (ice dynamics)
85      ! ========     
86      CALL lim_rhg
87      !
88      !
89      ! Control prints
90      IF(ln_ctl) THEN
91         CALL prt_ctl_info(' ')
92         CALL prt_ctl_info(' - Cell values : ')
93         CALL prt_ctl_info('   ~~~~~~~~~~~~~ ')
94         CALL prt_ctl(tab2d_1=divu_i    , clinfo1=' lim_dyn  : divu_i    :')
95         CALL prt_ctl(tab2d_1=delta_i   , clinfo1=' lim_dyn  : delta_i   :')
96         CALL prt_ctl(tab2d_1=strength  , clinfo1=' lim_dyn  : strength  :')
97         CALL prt_ctl(tab2d_1=e12t      , clinfo1=' lim_dyn  : cell area :')
98         CALL prt_ctl(tab2d_1=at_i      , clinfo1=' lim_dyn  : at_i      :')
99         CALL prt_ctl(tab2d_1=vt_i      , clinfo1=' lim_dyn  : vt_i      :')
100         CALL prt_ctl(tab2d_1=vt_s      , clinfo1=' lim_dyn  : vt_s      :')
101         CALL prt_ctl(tab2d_1=stress1_i , clinfo1=' lim_dyn  : stress1_i :')
102         CALL prt_ctl(tab2d_1=stress2_i , clinfo1=' lim_dyn  : stress2_i :')
103         CALL prt_ctl(tab2d_1=stress12_i, clinfo1=' lim_dyn  : stress12_i:')
104         DO jl = 1, jpl
105            CALL prt_ctl_info(' ')
106            CALL prt_ctl_info(' - Category : ', ivar1=jl)
107            CALL prt_ctl_info('   ~~~~~~~~~~')
108            CALL prt_ctl(tab2d_1=a_i   (:,:,jl)   , clinfo1= ' lim_dyn  : a_i      : ')
109            CALL prt_ctl(tab2d_1=ht_i  (:,:,jl)   , clinfo1= ' lim_dyn  : ht_i     : ')
110            CALL prt_ctl(tab2d_1=ht_s  (:,:,jl)   , clinfo1= ' lim_dyn  : ht_s     : ')
111            CALL prt_ctl(tab2d_1=v_i   (:,:,jl)   , clinfo1= ' lim_dyn  : v_i      : ')
112            CALL prt_ctl(tab2d_1=v_s   (:,:,jl)   , clinfo1= ' lim_dyn  : v_s      : ')
113            CALL prt_ctl(tab2d_1=e_s   (:,:,1,jl) , clinfo1= ' lim_dyn  : e_s      : ')
114            CALL prt_ctl(tab2d_1=t_su  (:,:,jl)   , clinfo1= ' lim_dyn  : t_su     : ')
115            CALL prt_ctl(tab2d_1=t_s   (:,:,1,jl) , clinfo1= ' lim_dyn  : t_snow   : ')
116            CALL prt_ctl(tab2d_1=sm_i  (:,:,jl)   , clinfo1= ' lim_dyn  : sm_i     : ')
117            CALL prt_ctl(tab2d_1=smv_i (:,:,jl)   , clinfo1= ' lim_dyn  : smv_i    : ')
118            DO jk = 1, nlay_i
119               CALL prt_ctl_info(' ')
120               CALL prt_ctl_info(' - Layer : ', ivar1=jk)
121               CALL prt_ctl_info('   ~~~~~~~')
122               CALL prt_ctl(tab2d_1=t_i(:,:,jk,jl) , clinfo1= ' lim_dyn  : t_i      : ')
123               CALL prt_ctl(tab2d_1=e_i(:,:,jk,jl) , clinfo1= ' lim_dyn  : e_i      : ')
124            END DO
125         END DO
126      ENDIF
127      !
128      ! conservation test
129      IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limdyn', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b)
130      !
131      IF( nn_timing == 1 )  CALL timing_stop('limdyn')
132
133   END SUBROUTINE lim_dyn
134
135
136   SUBROUTINE lim_dyn_init
137      !!-------------------------------------------------------------------
138      !!                  ***  ROUTINE lim_dyn_init  ***
139      !!
140      !! ** Purpose : Physical constants and parameters linked to the ice
141      !!      dynamics
142      !!
143      !! ** Method  :  Read the namicedyn namelist and check the ice-dynamic
144      !!       parameter values called at the first timestep (nit000)
145      !!
146      !! ** input   :   Namelist namicedyn
147      !!-------------------------------------------------------------------
148      INTEGER  ::   ios                 ! Local integer output status for namelist read
149      NAMELIST/namicedyn/ nn_icestr, ln_icestr_bvf, rn_pe_rdg, rn_pstar, rn_crhg, rn_cio,  rn_creepl, rn_ecc, &
150         &                nn_nevp, rn_relast, ln_landfast, rn_gamma, rn_icebfr
151      !!-------------------------------------------------------------------
152
153      REWIND( numnam_ice_ref )              ! Namelist namicedyn in reference namelist : Ice dynamics
154      READ  ( numnam_ice_ref, namicedyn, IOSTAT = ios, ERR = 901)
155901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicedyn in reference namelist', lwp )
156
157      REWIND( numnam_ice_cfg )              ! Namelist namicedyn in configuration namelist : Ice dynamics
158      READ  ( numnam_ice_cfg, namicedyn, IOSTAT = ios, ERR = 902 )
159902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicedyn in configuration namelist', lwp )
160      IF(lwm) WRITE ( numoni, namicedyn )
161     
162      IF(lwp) THEN                        ! control print
163         WRITE(numout,*)
164         WRITE(numout,*) 'lim_dyn_init : ice parameters for ice dynamics '
165         WRITE(numout,*) '~~~~~~~~~~~~'
166         WRITE(numout,*)'    ice strength parameterization (0=Hibler 1=Rothrock)         nn_icestr     = ', nn_icestr 
167         WRITE(numout,*)'    Including brine volume in ice strength comp.                ln_icestr_bvf = ', ln_icestr_bvf
168         WRITE(numout,*)'    Ratio of ridging work to PotEner change in ridging          rn_pe_rdg     = ', rn_pe_rdg 
169         WRITE(numout,*) '   drag coefficient for oceanic stress                         rn_cio        = ', rn_cio
170         WRITE(numout,*) '   first bulk-rheology parameter                               rn_pstar      = ', rn_pstar
171         WRITE(numout,*) '   second bulk-rhelogy parameter                               rn_crhg       = ', rn_crhg
172         WRITE(numout,*) '   creep limit                                                 rn_creepl     = ', rn_creepl
173         WRITE(numout,*) '   eccentricity of the elliptical yield curve                  rn_ecc        = ', rn_ecc
174         WRITE(numout,*) '   number of iterations for subcycling                         nn_nevp       = ', nn_nevp
175         WRITE(numout,*) '   ratio of elastic timescale over ice time step               rn_relast     = ', rn_relast
176         WRITE(numout,*) '   Landfast: param (T or F)                                    ln_landfast   = ', ln_landfast
177         WRITE(numout,*) '   Landfast: fraction of ocean depth that ice must reach       rn_gamma      = ', rn_gamma
178         WRITE(numout,*) '   Landfast: maximum bottom stress per unit area of contact    rn_icebfr     = ', rn_icebfr
179      ENDIF
180      !
181      usecc2 = 1._wp / ( rn_ecc * rn_ecc )
182      rhoco  = rau0  * rn_cio
183      !
184   END SUBROUTINE lim_dyn_init
185
186#else
187   !!----------------------------------------------------------------------
188   !!   Default option          Empty module           NO LIM sea-ice model
189   !!----------------------------------------------------------------------
190CONTAINS
191   SUBROUTINE lim_dyn         ! Empty routine
192   END SUBROUTINE lim_dyn
193#endif 
194
195   !!======================================================================
196END MODULE limdyn
Note: See TracBrowser for help on using the repository browser.