Index: NEMO/branches/2018/dev_r9947_SI3_advection/cfgs/SHARED/field_def_nemo-ice.xml
===================================================================
--- NEMO/branches/2018/dev_r9947_SI3_advection/cfgs/SHARED/field_def_nemo-ice.xml (revision 10281)
+++ NEMO/branches/2018/dev_r9947_SI3_advection/cfgs/SHARED/field_def_nemo-ice.xml (revision 10312)
@@ -29,5 +29,6 @@
-
+
+
Index: NEMO/branches/2018/dev_r9947_SI3_advection/cfgs/SHARED/namelist_ice_ref
===================================================================
--- NEMO/branches/2018/dev_r9947_SI3_advection/cfgs/SHARED/namelist_ice_ref (revision 10281)
+++ NEMO/branches/2018/dev_r9947_SI3_advection/cfgs/SHARED/namelist_ice_ref (revision 10312)
@@ -52,14 +52,16 @@
ln_dynRHGADV = .false. ! dyn.: no ridge/raft & no corrections (rheology + advection)
ln_dynADV = .false. ! dyn.: only advection w prescribed vel.(rn_uvice + advection)
- rn_uice = 0.00001 ! prescribed ice u-velocity
- rn_vice = 0. ! prescribed ice v-velocity
- rn_ishlat = 2. ! free slip (0) ; partial slip (0-2) ; no slip (2) ; strong slip (>2)
- ln_landfast = .false. ! landfast ice parameterization (T or F)
- rn_gamma = 0.15 ! fraction of ocean depth that ice must reach to initiate landfast
- ! recommended range: [0.1 ; 0.25]
- rn_icebfr = 10. ! maximum bottom stress per unit area of contact [N/m2]
- ! a very large value ensures ice velocity=0 even with a small contact area
- ! recommended range: ?? (should be greater than atm-ice stress => >0.1 N/m2)
- rn_lfrelax = 1.e-5 ! relaxation time scale to reach static friction [s-1]
+ rn_uice = 0.5 ! prescribed ice u-velocity
+ rn_vice = 0.5 ! prescribed ice v-velocity
+ rn_ishlat = 2. ! lbc : free slip (0) ; partial slip (0-2) ; no slip (2) ; strong slip (>2)
+ ln_landfast_L16 = .false. ! landfast: parameterization from Lemieux 2016
+ ln_landfast_home = .false. ! landfast: parameterization from "home made"
+ rn_depfra = 0.125 ! fraction of ocean depth that ice must reach to initiate landfast
+ ! recommended range: [0.1 ; 0.25] - L16=0.125 - home=0.15
+ rn_icebfr = 15. ! ln_landfast_L16: maximum bottom stress per unit volume [N/m3]
+ ! ln_landfast_home: maximum bottom stress per unit area of contact [N/m2]
+ ! recommended range: ?? L16=15 - home=10
+ rn_lfrelax = 1.e-5 ! relaxation time scale to reach static friction [s-1]
+ rn_tensile = 0.2 ! isotropic tensile strength
/
!------------------------------------------------------------------------------
Index: NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/ice.F90
===================================================================
--- NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/ice.F90 (revision 10281)
+++ NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/ice.F90 (revision 10312)
@@ -129,10 +129,15 @@
! !!** ice-dynamics namelist (namdyn) **
REAL(wp), PUBLIC :: rn_ishlat !: lateral boundary condition for sea-ice
- LOGICAL , PUBLIC :: ln_landfast !: landfast ice parameterization (T or F)
- REAL(wp), PUBLIC :: rn_gamma !: fraction of ocean depth that ice must reach to initiate landfast ice
- REAL(wp), PUBLIC :: rn_icebfr !: maximum bottom stress per unit area of contact (landfast ice)
- REAL(wp), PUBLIC :: rn_lfrelax !: relaxation time scale (s-1) to reach static friction (landfast ice)
- !
- ! !!** ice-rheology namelist (namrhg) **
+ LOGICAL , PUBLIC :: ln_landfast_L16 !: landfast ice parameterizationfrom lemieux2016
+ LOGICAL , PUBLIC :: ln_landfast_home !: landfast ice parameterizationfrom home made
+ REAL(wp), PUBLIC :: rn_depfra !: fraction of ocean depth that ice must reach to initiate landfast ice
+ REAL(wp), PUBLIC :: rn_icebfr !: maximum bottom stress per unit area of contact (lemieux2016) or per unit volume (home)
+ REAL(wp), PUBLIC :: rn_lfrelax !: relaxation time scale (s-1) to reach static friction
+ REAL(wp), PUBLIC :: rn_tensile !: isotropic tensile strength
+ !
+ ! !!** ice-ridging/rafting namelist (namdyn_rdgrft) **
+ REAL(wp), PUBLIC :: rn_crhg !: determines changes in ice strength (also used for landfast param)
+ !
+ ! !!** ice-rheology namelist (namdyn_rhg) **
LOGICAL , PUBLIC :: ln_aEVP !: using adaptive EVP (T or F)
REAL(wp), PUBLIC :: rn_creepl !: creep limit : has to be under 1.0e-9
Index: NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn.F90
===================================================================
--- NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn.F90 (revision 10281)
+++ NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn.F90 (revision 10312)
@@ -84,14 +84,12 @@
WRITE(numout,*)'~~~~~~~'
ENDIF
-
!
- IF( ln_landfast ) THEN !-- Landfast ice parameterization: define max bottom friction
+ IF( ln_landfast_home ) THEN !-- Landfast ice parameterization
tau_icebfr(:,:) = 0._wp
DO jl = 1, jpl
- WHERE( h_i(:,:,jl) > ht_n(:,:) * rn_gamma ) tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr
+ WHERE( h_i(:,:,jl) > ht_n(:,:) * rn_depfra ) tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr
END DO
- IF( iom_use('tau_icebfr') ) CALL iom_put( 'tau_icebfr', tau_icebfr )
- ENDIF
-
+ ENDIF
+ !
zhmax(:,:,:) = h_i_b(:,:,:) !-- Record max of the surrounding 9-pts ice thick. (for CALL Hbig)
DO jl = 1, jpl
@@ -244,6 +242,6 @@
!!
NAMELIST/namdyn/ ln_dynFULL, ln_dynRHGADV, ln_dynADV, rn_uice, rn_vice, &
- & rn_ishlat , &
- & ln_landfast, rn_gamma , rn_icebfr, rn_lfrelax
+ & rn_ishlat , &
+ & ln_landfast_L16, ln_landfast_home, rn_depfra, rn_icebfr, rn_lfrelax, rn_tensile
!!-------------------------------------------------------------------
!
@@ -261,13 +259,15 @@
WRITE(numout,*) '~~~~~~~~~~~~'
WRITE(numout,*) ' Namelist namdyn:'
- WRITE(numout,*) ' Full ice dynamics (rhg + adv + ridge/raft + corr) ln_dynFULL = ', ln_dynFULL
- WRITE(numout,*) ' No ridge/raft & No cor (rhg + adv) ln_dynRHGADV = ', ln_dynRHGADV
- WRITE(numout,*) ' Advection only (rn_uvice + adv) ln_dynADV = ', ln_dynADV
- WRITE(numout,*) ' with prescribed velocity given by (u,v)_ice = (rn_uice,rn_vice) = (', rn_uice,',', rn_vice,')'
- WRITE(numout,*) ' lateral boundary condition for sea ice dynamics rn_ishlat = ', rn_ishlat
- WRITE(numout,*) ' Landfast: param (T or F) ln_landfast = ', ln_landfast
- WRITE(numout,*) ' fraction of ocean depth that ice must reach rn_gamma = ', rn_gamma
- WRITE(numout,*) ' maximum bottom stress per unit area of contact rn_icebfr = ', rn_icebfr
- WRITE(numout,*) ' relax time scale (s-1) to reach static friction rn_lfrelax = ', rn_lfrelax
+ WRITE(numout,*) ' Full ice dynamics (rhg + adv + ridge/raft + corr) ln_dynFULL = ', ln_dynFULL
+ WRITE(numout,*) ' No ridge/raft & No cor (rhg + adv) ln_dynRHGADV = ', ln_dynRHGADV
+ WRITE(numout,*) ' Advection only (rn_uvice + adv) ln_dynADV = ', ln_dynADV
+ WRITE(numout,*) ' with prescribed velocity given by (u,v)_ice = (rn_uice,rn_vice) = (', rn_uice,',', rn_vice,')'
+ WRITE(numout,*) ' lateral boundary condition for sea ice dynamics rn_ishlat = ', rn_ishlat
+ WRITE(numout,*) ' Landfast: param from Lemieux 2016 ln_landfast_L16 = ', ln_landfast_L16
+ WRITE(numout,*) ' Landfast: param from home made ln_landfast_home= ', ln_landfast_home
+ WRITE(numout,*) ' fraction of ocean depth that ice must reach rn_depfra = ', rn_depfra
+ WRITE(numout,*) ' maximum bottom stress per unit area of contact rn_icebfr = ', rn_icebfr
+ WRITE(numout,*) ' relax time scale (s-1) to reach static friction rn_lfrelax = ', rn_lfrelax
+ WRITE(numout,*) ' isotropic tensile strength rn_tensile = ', rn_tensile
WRITE(numout,*)
ENDIF
@@ -289,6 +289,10 @@
ELSEIF ( 2. < rn_ishlat ) THEN ; IF(lwp) WRITE(numout,*) ' ===>>> ice lateral strong-slip'
ENDIF
- ! !--- NO Landfast ice : set to zero once for all
- IF( .NOT.ln_landfast ) tau_icebfr(:,:) = 0._wp
+ ! !--- Landfast ice
+ IF( .NOT.ln_landfast_L16 .AND. .NOT.ln_landfast_home ) tau_icebfr(:,:) = 0._wp
+ !
+ IF ( ln_landfast_L16 .AND. ln_landfast_home ) THEN
+ CALL ctl_stop( 'ice_dyn_init: choose one and only one landfast parameterization (ln_landfast_L16 or ln_landfast_home)' )
+ ENDIF
!
CALL ice_dyn_rdgrft_init ! set ice ridging/rafting parameters
Index: NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_rdgrft.F90
===================================================================
--- NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_rdgrft.F90 (revision 10281)
+++ NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_rdgrft.F90 (revision 10312)
@@ -39,15 +39,15 @@
! Variables shared among ridging subroutines
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: closing_net ! net rate at which area is removed (1/s)
- ! ! (ridging ice area - area of new ridges) / dt
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: opning ! rate of opening due to divergence/shear
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: closing_gross ! rate at which area removed, not counting area of new ridges
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: apartf ! participation function; fraction of ridging/closing associated w/ category n
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hrmin ! minimum ridge thickness
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hrmax ! maximum ridge thickness
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hraft ! thickness of rafted ice
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hi_hrdg ! thickness of ridging ice / mean ridge thickness
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: aridge ! participating ice ridging
- REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: araft ! participating ice rafting
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: closing_net ! net rate at which area is removed (1/s)
+ ! ! (ridging ice area - area of new ridges) / dt
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: opning ! rate of opening due to divergence/shear
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) :: closing_gross ! rate at which area removed, not counting area of new ridges
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: apartf ! participation function; fraction of ridging/closing associated w/ category n
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hrmin ! minimum ridge thickness
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hrmax ! maximum ridge thickness
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hraft ! thickness of rafted ice
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: hi_hrdg ! thickness of ridging ice / mean ridge thickness
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: aridge ! participating ice ridging
+ REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) :: araft ! participating ice rafting
REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ze_i_2d
REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: ze_s_2d
@@ -59,5 +59,4 @@
LOGICAL :: ln_str_H79 ! ice strength parameterization (Hibler79)
REAL(wp) :: rn_pstar ! determines ice strength, Hibler JPO79
- REAL(wp) :: rn_crhg ! determines changes in ice strength
REAL(wp) :: rn_csrdg ! fraction of shearing energy contributing to ridging
LOGICAL :: ln_partf_lin ! participation function linear (Thorndike et al. (1975))
Index: NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_rhg_evp.F90
===================================================================
--- NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_rhg_evp.F90 (revision 10281)
+++ NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_rhg_evp.F90 (revision 10312)
@@ -118,7 +118,9 @@
REAL(wp) :: ecc2, z1_ecc2 ! square of yield ellipse eccenticity
REAL(wp) :: zalph1, z1_alph1, zalph2, z1_alph2 ! alpha coef from Bouillon 2009 or Kimmritz 2017
- REAL(wp) :: zm1, zm2, zm3, zmassU, zmassV ! ice/snow mass
+ REAL(wp) :: zm1, zm2, zm3, zmassU, zmassV, zvU, zvV ! ice/snow mass and volume
REAL(wp) :: zdelta, zp_delf, zds2, zdt, zdt2, zdiv, zdiv2 ! temporary scalars
REAL(wp) :: zTauO, zTauB, zTauE, zvel ! temporary scalars
+ REAL(wp) :: zkt ! isotropic tensile strength for landfast ice
+ REAL(wp) :: zvCr ! critical ice volume above which ice is landfast
!
REAL(wp) :: zresm ! Maximal error on ice velocity
@@ -136,4 +138,5 @@
REAL(wp), DIMENSION(jpi,jpj) :: zmf ! coriolis parameter at T points
REAL(wp), DIMENSION(jpi,jpj) :: zTauU_ia , ztauV_ia ! ice-atm. stress at U-V points
+ REAL(wp), DIMENSION(jpi,jpj) :: zTauU_ib , ztauV_ib ! ice-bottom stress at U-V points (landfast param)
REAL(wp), DIMENSION(jpi,jpj) :: zspgU , zspgV ! surface pressure gradient at U/V points
REAL(wp), DIMENSION(jpi,jpj) :: v_oceU, u_oceV, v_iceU, u_iceV ! ocean/ice u/v component on V/U points
@@ -255,5 +258,9 @@
END DO
END DO
-
+
+ ! landfast param from Lemieux(2016): add isotropic tensile strength (following Konig Beatty and Holland, 2010)
+ IF( ln_landfast_L16 .OR. ln_landfast_home ) THEN ; zkt = rn_tensile
+ ELSE ; zkt = 0._wp
+ ENDIF
!
!------------------------------------------------------------------------------!
@@ -273,5 +280,5 @@
zpice(:,:) = ssh_m(:,:) + ( zintn * snwice_mass(:,:) + zintb * snwice_mass_b(:,:) ) * r1_rau0
!
- ELSE !== non-embedded sea ice: use ocean surface for slope calculation ==!
+ ELSE !== non-embedded sea ice: use ocean surface for slope calculation ==!
zpice(:,:) = ssh_m(:,:)
ENDIF
@@ -328,4 +335,43 @@
CALL lbc_lnk_multi( zmf, 'T', 1., zdt_m, 'T', 1. )
!
+ ! !== Landfast ice parameterization ==!
+ !
+ IF( ln_landfast_L16 ) THEN !-- Lemieux 2016
+ DO jj = 2, jpjm1
+ DO ji = fs_2, fs_jpim1
+ ! ice thickness at U-V points
+ zvU = 0.5_wp * ( vt_i(ji,jj) * e1e2t(ji,jj) + vt_i(ji+1,jj) * e1e2t(ji+1,jj) ) * r1_e1e2u(ji,jj) * umask(ji,jj,1)
+ zvV = 0.5_wp * ( vt_i(ji,jj) * e1e2t(ji,jj) + vt_i(ji,jj+1) * e1e2t(ji,jj+1) ) * r1_e1e2v(ji,jj) * vmask(ji,jj,1)
+ ! ice-bottom stress at U points
+ zvCr = zaU(ji,jj) * rn_depfra * hu_n(ji,jj)
+ zTauU_ib(ji,jj) = rn_icebfr * MAX( 0._wp, zvU - zvCr ) * EXP( -rn_crhg * ( 1._wp - zaU(ji,jj) ) )
+ ! ice-bottom stress at V points
+ zvCr = zaV(ji,jj) * rn_depfra * hv_n(ji,jj)
+ zTauV_ib(ji,jj) = rn_icebfr * MAX( 0._wp, zvV - zvCr ) * EXP( -rn_crhg * ( 1._wp - zaV(ji,jj) ) )
+ ! ice_bottom stress at T points
+ zvCr = at_i(ji,jj) * rn_depfra * ht_n(ji,jj)
+ tau_icebfr(ji,jj) = rn_icebfr * MAX( 0._wp, vt_i(ji,jj) - zvCr ) * EXP( -rn_crhg * ( 1._wp - at_i(ji,jj) ) )
+ END DO
+ END DO
+ CALL lbc_lnk( tau_icebfr(:,:), 'T', 1. )
+ !
+ ELSEIF( ln_landfast_home ) THEN !-- Home made
+ DO jj = 2, jpjm1
+ DO ji = fs_2, fs_jpim1
+ zTauU_ib(ji,jj) = tau_icebfr(ji,jj)
+ zTauV_ib(ji,jj) = tau_icebfr(ji,jj)
+ END DO
+ END DO
+ !
+ ELSE !-- no landfast
+ DO jj = 2, jpjm1
+ DO ji = fs_2, fs_jpim1
+ zTauU_ib(ji,jj) = 0._wp
+ zTauV_ib(ji,jj) = 0._wp
+ END DO
+ END DO
+ ENDIF
+ IF( iom_use('tau_icebfr') ) CALL iom_put( 'tau_icebfr', tau_icebfr(:,:) )
+
!------------------------------------------------------------------------------!
! 3) Solution of the momentum equation, iterative procedure
@@ -391,7 +437,7 @@
ENDIF
- ! stress at T points
- zs1(ji,jj) = ( zs1(ji,jj) * zalph1 + zp_delt(ji,jj) * ( zdiv - zdelta ) ) * z1_alph1
- zs2(ji,jj) = ( zs2(ji,jj) * zalph2 + zp_delt(ji,jj) * ( zdt * z1_ecc2 ) ) * z1_alph2
+ ! stress at T points (zkt/=0 if landfast)
+ zs1(ji,jj) = ( zs1(ji,jj) * zalph1 + zp_delt(ji,jj) * ( zdiv * (1._wp + zkt) - zdelta * (1._wp - zkt) ) ) * z1_alph1
+ zs2(ji,jj) = ( zs2(ji,jj) * zalph2 + zp_delt(ji,jj) * ( zdt * z1_ecc2 * (1._wp + zkt) ) ) * z1_alph2
END DO
@@ -412,6 +458,6 @@
zp_delf = 0.25_wp * ( zp_delt(ji,jj) + zp_delt(ji+1,jj) + zp_delt(ji,jj+1) + zp_delt(ji+1,jj+1) )
- ! stress at F points
- zs12(ji,jj)= ( zs12(ji,jj) * zalph2 + zp_delf * ( zds(ji,jj) * z1_ecc2 ) * 0.5_wp ) * z1_alph2
+ ! stress at F points (zkt/=0 if landfast)
+ zs12(ji,jj)= ( zs12(ji,jj) * zalph2 + zp_delf * ( zds(ji,jj) * z1_ecc2 * (1._wp + zkt) ) * 0.5_wp ) * z1_alph2
END DO
@@ -461,6 +507,6 @@
!
! !--- tau_bottom/v_ice
- zvel = MAX( zepsi, SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) )
- zTauB = - tau_icebfr(ji,jj) / zvel
+ zvel = 5.e-05_wp + SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) )
+ zTauB = - zTauV_ib(ji,jj) / zvel
!
! !--- Coriolis at V-points (energy conserving formulation)
@@ -473,5 +519,5 @@
!
! !--- landfast switch => 0 = static friction ; 1 = sliding friction
- rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - tau_icebfr(ji,jj) ) - SIGN( 1._wp, zTauE ) ) )
+ rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - zTauV_ib(ji,jj) ) - SIGN( 1._wp, zTauE ) ) )
!
IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017)
@@ -509,6 +555,6 @@
!
! !--- tau_bottom/u_ice
- zvel = MAX( zepsi, SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) )
- zTauB = - tau_icebfr(ji,jj) / zvel
+ zvel = 5.e-05_wp + SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) )
+ zTauB = - zTauU_ib(ji,jj) / zvel
!
! !--- Coriolis at U-points (energy conserving formulation)
@@ -521,5 +567,5 @@
!
! !--- landfast switch => 0 = static friction ; 1 = sliding friction
- rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - tau_icebfr(ji,jj) ) - SIGN( 1._wp, zTauE ) ) )
+ rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - zTauU_ib(ji,jj) ) - SIGN( 1._wp, zTauE ) ) )
!
IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017)
@@ -559,6 +605,6 @@
!
! !--- tau_bottom/u_ice
- zvel = MAX( zepsi, SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) )
- zTauB = - tau_icebfr(ji,jj) / zvel
+ zvel = 5.e-05_wp + SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) )
+ zTauB = - zTauU_ib(ji,jj) / zvel
!
! !--- Coriolis at U-points (energy conserving formulation)
@@ -571,5 +617,5 @@
!
! !--- landfast switch => 0 = static friction ; 1 = sliding friction
- rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - tau_icebfr(ji,jj) ) - SIGN( 1._wp, zTauE ) ) )
+ rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - zTauU_ib(ji,jj) ) - SIGN( 1._wp, zTauE ) ) )
!
IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017)
@@ -607,6 +653,6 @@
!
! !--- tau_bottom/v_ice
- zvel = MAX( zepsi, SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) )
- ztauB = - tau_icebfr(ji,jj) / zvel
+ zvel = 5.e-05_wp + SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) )
+ zTauB = - zTauV_ib(ji,jj) / zvel
!
! !--- Coriolis at v-points (energy conserving formulation)
@@ -619,5 +665,5 @@
!
! !--- landfast switch => 0 = static friction ; 1 = sliding friction
- rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, zTauE - tau_icebfr(ji,jj) ) - SIGN( 1._wp, zTauE ) ) )
+ rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, zTauE - zTauV_ib(ji,jj) ) - SIGN( 1._wp, zTauE ) ) )
!
IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017)
Index: NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icewri.F90
===================================================================
--- NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icewri.F90 (revision 10281)
+++ NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icewri.F90 (revision 10312)
@@ -50,5 +50,5 @@
INTEGER :: ji, jj, jk, jl ! dummy loop indices
REAL(wp) :: z2da, z2db, zrho1, zrho2
- REAL(wp), DIMENSION(jpi,jpj) :: z2d ! 2D workspace
+ REAL(wp), DIMENSION(jpi,jpj) :: z2d, zfast ! 2D workspace
REAL(wp), DIMENSION(jpi,jpj) :: zmsk00, zmsk05, zmsk15, zmsksn ! O%, 5% and 15% concentration mask and snow mask
REAL(wp), DIMENSION(jpi,jpj,jpl) :: zmsk00l, zmsksnl ! cat masks
@@ -132,5 +132,6 @@
IF( iom_use('vtau_ai' ) ) CALL iom_put( "vtau_ai", vtau_ice * zmsk00 ) ! Wind stress term in force balance (y)
- IF( iom_use('icevel') ) THEN
+ IF( iom_use('icevel') .OR. iom_use('fasticepres') ) THEN
+ ! module of ice velocity
DO jj = 2 , jpjm1
DO ji = 2 , jpim1
@@ -141,5 +142,11 @@
END DO
CALL lbc_lnk( z2d, 'T', 1. )
- IF( iom_use('icevel') ) CALL iom_put( "icevel" , z2d ) ! ice velocity module
+ IF( iom_use('icevel') ) CALL iom_put( "icevel" , z2d )
+
+ ! record presence of fast ice
+ WHERE( z2d(:,:) < 5.e-04_wp .AND. zmsk00(:,:) == 1._wp ) ; zfast(:,:) = 1._wp
+ ELSEWHERE ; zfast(:,:) = 0._wp
+ END WHERE
+ IF( iom_use('fasticepres') ) CALL iom_put( "fasticepres" , zfast )
ENDIF