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.
icerhg.F90 in branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3 – NEMO

source: branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icerhg.F90 @ 8407

Last change on this file since 8407 was 8407, checked in by clem, 7 years ago

change calls in icestp.F90 for rheology
File size: 8.6 KB
Line 
1MODULE icerhg
2   !!======================================================================
3   !!                     ***  MODULE  icerhg  ***
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   !!    ice_rhg      : computes ice velocities
15   !!    ice_rhg_init : initialization and namelist read
16   !!----------------------------------------------------------------------
17   USE phycst           ! physical constants
18   USE dom_oce          ! ocean space and time domain
19   USE ice              ! LIM-3 variables
20   USE icerhg_evp       ! EVP rheology
21   USE limcons          ! conservation tests
22   USE limctl           ! control prints
23   USE limvar
24   !
25   USE lbclnk           ! lateral boundary conditions - MPP exchanges
26   USE lib_mpp          ! MPP library
27   USE in_out_manager   ! I/O manager
28   USE lib_fortran      ! glob_sum
29   USE timing           ! Timing
30
31   IMPLICIT NONE
32   PRIVATE
33
34   PUBLIC   ice_rhg        ! routine called by icestp.F90
35   PUBLIC   ice_rhg_init   ! routine called by icestp.F90
36
37   !! * Substitutions
38#  include "vectopt_loop_substitute.h90"
39   !!----------------------------------------------------------------------
40   !! NEMO/LIM3 4.0 , UCL - NEMO Consortium (2011)
41   !! $Id: icerhg.F90 8378 2017-07-26 13:55:59Z clem $
42   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt)
43   !!----------------------------------------------------------------------
44CONTAINS
45
46   SUBROUTINE ice_rhg( kt )
47      !!-------------------------------------------------------------------
48      !!               ***  ROUTINE ice_rhg  ***
49      !!               
50      !! ** Purpose :   compute ice velocity
51      !!
52      !! ** Action  : comupte - ice velocity (u_ice, v_ice)
53      !!                      - 3 components of the stress tensor (stress1_i, stress2_i, stress12_i)
54      !!                      - shear, divergence and delta (shear_i, divu_i, delta_i)
55      !!--------------------------------------------------------------------
56      INTEGER, INTENT(in) ::   kt     ! number of iteration
57      !!
58      INTEGER  :: jl ! dummy loop indices
59      REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 
60      !!--------------------------------------------------------------------
61
62      IF( nn_timing == 1 )  CALL timing_start('icerhg')
63
64      CALL lim_var_agg(1)   ! aggregate ice categories
65      !
66      ! conservation test
67      IF( ln_limdiachk ) CALL lim_cons_hsm(0, 'icerhg', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b)
68     
69      ! Landfast ice parameterization: define max bottom friction
70      tau_icebfr(:,:) = 0._wp
71      IF( ln_landfast ) THEN
72         DO jl = 1, jpl
73            WHERE( ht_i(:,:,jl) > ht_n(:,:) * rn_gamma )  tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr
74         END DO
75      ENDIF
76     
77      ! -----------------------
78      ! Rheology (ice dynamics)
79      ! -----------------------   
80      IF( nn_limdyn /= 0 ) THEN                          ! -- Ice dynamics
81
82         CALL ice_rhg_evp( stress1_i, stress2_i, stress12_i, u_ice, v_ice, shear_i, divu_i, delta_i )
83
84      ELSE
85
86         u_ice(:,:) = rn_uice * umask(:,:,1)             !     or prescribed velocity
87         v_ice(:,:) = rn_vice * vmask(:,:,1)
88         !!CALL RANDOM_NUMBER(u_ice(:,:))
89         !!CALL RANDOM_NUMBER(v_ice(:,:))
90
91      ENDIF
92      !
93      ! conservation test
94      IF( ln_limdiachk ) CALL lim_cons_hsm(1, 'icerhg', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b)
95
96      ! Control prints
97      IF( ln_ctl )       CALL lim_prt3D( 'icerhg' )
98      !
99      IF( nn_timing == 1 )  CALL timing_stop('icerhg')
100
101   END SUBROUTINE ice_rhg
102
103
104   SUBROUTINE ice_rhg_init
105      !!-------------------------------------------------------------------
106      !!                  ***  ROUTINE ice_rhg_init  ***
107      !!
108      !! ** Purpose : Physical constants and parameters linked to the ice
109      !!      dynamics
110      !!
111      !! ** Method  :  Read the namicedyn namelist and check the ice-dynamic
112      !!       parameter values called at the first timestep (nit000)
113      !!
114      !! ** input   :   Namelist namicedyn
115      !!-------------------------------------------------------------------
116      INTEGER  ::   ios                 ! Local integer output status for namelist read
117      NAMELIST/namicedyn/ nn_limadv, nn_limadv_ord,                                &
118         &                nn_icestr, rn_pe_rdg, rn_pstar, rn_crhg, ln_icestr_bvf,  &
119         &                rn_ishlat, rn_cio, rn_creepl, rn_ecc,                    &
120         &                nn_nevp, rn_relast, ln_landfast, rn_gamma, rn_icebfr, rn_lfrelax
121      !!-------------------------------------------------------------------
122
123      REWIND( numnam_ice_ref )              ! Namelist namicedyn in reference namelist : Ice dynamics
124      READ  ( numnam_ice_ref, namicedyn, IOSTAT = ios, ERR = 901)
125901   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicedyn in reference namelist', lwp )
126
127      REWIND( numnam_ice_cfg )              ! Namelist namicedyn in configuration namelist : Ice dynamics
128      READ  ( numnam_ice_cfg, namicedyn, IOSTAT = ios, ERR = 902 )
129902   IF( ios /= 0 ) CALL ctl_nam ( ios , 'namicedyn in configuration namelist', lwp )
130      IF(lwm) WRITE ( numoni, namicedyn )
131     
132      IF(lwp) THEN                        ! control print
133         WRITE(numout,*)
134         WRITE(numout,*) 'ice_rhg_init : ice parameters for ice dynamics '
135         WRITE(numout,*) '~~~~~~~~~~~~'
136         ! limtrp
137         WRITE(numout,*)'    choose the advection scheme (-1=Prather, 0=Ulimate-Macho)   nn_limadv     = ', nn_limadv 
138         WRITE(numout,*)'    choose the order of the scheme (if ultimate)                nn_limadv_ord = ', nn_limadv_ord 
139         ! limitd_me
140         WRITE(numout,*)'    ice strength parameterization (0=Hibler 1=Rothrock)         nn_icestr     = ', nn_icestr 
141         WRITE(numout,*)'    Ratio of ridging work to PotEner change in ridging          rn_pe_rdg     = ', rn_pe_rdg 
142         WRITE(numout,*) '   first bulk-rheology parameter                               rn_pstar      = ', rn_pstar
143         WRITE(numout,*) '   second bulk-rhelogy parameter                               rn_crhg       = ', rn_crhg
144         WRITE(numout,*)'    Including brine volume in ice strength comp.                ln_icestr_bvf = ', ln_icestr_bvf
145         ! icerhg_evp
146         WRITE(numout,*) '   lateral boundary condition for sea ice dynamics             rn_ishlat     = ', rn_ishlat
147         WRITE(numout,*) '   drag coefficient for oceanic stress                         rn_cio        = ', rn_cio
148         WRITE(numout,*) '   creep limit                                                 rn_creepl     = ', rn_creepl
149         WRITE(numout,*) '   eccentricity of the elliptical yield curve                  rn_ecc        = ', rn_ecc
150         WRITE(numout,*) '   number of iterations for subcycling                         nn_nevp       = ', nn_nevp
151         WRITE(numout,*) '   ratio of elastic timescale over ice time step               rn_relast     = ', rn_relast
152         WRITE(numout,*) '   Landfast: param (T or F)                                    ln_landfast   = ', ln_landfast
153         WRITE(numout,*) '      T: fraction of ocean depth that ice must reach           rn_gamma      = ', rn_gamma
154         WRITE(numout,*) '      T: maximum bottom stress per unit area of contact        rn_icebfr     = ', rn_icebfr
155         WRITE(numout,*) '      T: relax time scale (s-1) to reach static friction       rn_lfrelax    = ', rn_lfrelax
156      ENDIF
157      !
158      IF     (      rn_ishlat == 0.                ) THEN   ;   IF(lwp) WRITE(numout,*) '   ice lateral  free-slip '
159      ELSEIF (      rn_ishlat == 2.                ) THEN   ;   IF(lwp) WRITE(numout,*) '   ice lateral  no-slip '
160      ELSEIF ( 0. < rn_ishlat .AND. rn_ishlat < 2. ) THEN   ;   IF(lwp) WRITE(numout,*) '   ice lateral  partial-slip '
161      ELSEIF ( 2. < rn_ishlat                      ) THEN   ;   IF(lwp) WRITE(numout,*) '   ice lateral  strong-slip '
162      ENDIF
163      !
164   END SUBROUTINE ice_rhg_init
165
166#endif 
167
168   !!======================================================================
169END MODULE icerhg
Note: See TracBrowser for help on using the repository browser.