Changeset 10312


Ignore:
Timestamp:
2018-11-15T10:46:57+01:00 (23 months ago)
Author:
clem
Message:

add the parameterization of landfast ice from Lemieux2015 and 2016 with the addition of isotropic tensile strength

Location:
NEMO/branches/2018/dev_r9947_SI3_advection
Files:
7 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2018/dev_r9947_SI3_advection/cfgs/SHARED/field_def_nemo-ice.xml

    r9943 r10312  
    2929          <field id="icemask15"    long_name="Ice mask (0 if ice conc. lower than 15%, 1 otherwise)"   standard_name="sea_ice_mask15"                            unit="" /> 
    3030     <field id="icepres"      long_name="Fraction of time steps with sea ice"                     standard_name="sea_ice_time_fraction"                     unit="" /> 
    31       
     31          <field id="fasticepres"  long_name="Fraction of time steps with landfast ice"                standard_name="fast_ice_time_fraction"                    unit="" /> 
     32  
    3233     <!-- general fields --> 
    3334          <field id="icemass"      long_name="Sea-ice mass per area"                                   standard_name="sea_ice_amount"                            unit="kg/m2"/> 
  • NEMO/branches/2018/dev_r9947_SI3_advection/cfgs/SHARED/namelist_ice_ref

    r9801 r10312  
    5252   ln_dynRHGADV     = .false.         !  dyn.: no ridge/raft & no corrections  (rheology + advection) 
    5353   ln_dynADV        = .false.         !  dyn.: only advection w prescribed vel.(rn_uvice + advection) 
    54       rn_uice       =   0.00001       !        prescribed ice u-velocity 
    55       rn_vice       =   0.            !        prescribed ice v-velocity 
    56    rn_ishlat        =   2.            !  free slip (0) ; partial slip (0-2) ; no slip (2) ; strong slip (>2) 
    57    ln_landfast      = .false.         !  landfast ice parameterization (T or F)                            
    58       rn_gamma      =   0.15          !     fraction of ocean depth that ice must reach to initiate landfast 
    59                                       !        recommended range: [0.1 ; 0.25] 
    60       rn_icebfr     =  10.            !     maximum bottom stress per unit area of contact [N/m2]                  
    61                                       !        a very large value ensures ice velocity=0 even with a small contact area 
    62                                       !        recommended range: ?? (should be greater than atm-ice stress => >0.1 N/m2) 
    63       rn_lfrelax    =   1.e-5         !     relaxation time scale to reach static friction [s-1] 
     54      rn_uice       =   0.5           !        prescribed ice u-velocity 
     55      rn_vice       =   0.5           !        prescribed ice v-velocity 
     56   rn_ishlat        =   2.            !  lbc : free slip (0) ; partial slip (0-2) ; no slip (2) ; strong slip (>2) 
     57   ln_landfast_L16  = .false.         !  landfast: parameterization from Lemieux 2016 
     58   ln_landfast_home = .false.         !  landfast: parameterization from "home made" 
     59      rn_depfra     =   0.125         !        fraction of ocean depth that ice must reach to initiate landfast 
     60                                      !          recommended range: [0.1 ; 0.25] - L16=0.125 - home=0.15 
     61      rn_icebfr     =  15.            !        ln_landfast_L16:  maximum bottom stress per unit volume [N/m3] 
     62                                      !        ln_landfast_home: maximum bottom stress per unit area of contact [N/m2] 
     63                                      !          recommended range: ?? L16=15 - home=10 
     64      rn_lfrelax    =   1.e-5         !        relaxation time scale to reach static friction [s-1] 
     65      rn_tensile    =   0.2           !        isotropic tensile strength 
    6466/ 
    6567!------------------------------------------------------------------------------ 
  • NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/ice.F90

    r9922 r10312  
    129129   !                                     !!** ice-dynamics namelist (namdyn) ** 
    130130   REAL(wp), PUBLIC ::   rn_ishlat        !: lateral boundary condition for sea-ice 
    131    LOGICAL , PUBLIC ::   ln_landfast      !: landfast ice parameterization (T or F)  
    132    REAL(wp), PUBLIC ::   rn_gamma         !:    fraction of ocean depth that ice must reach to initiate landfast ice 
    133    REAL(wp), PUBLIC ::   rn_icebfr        !:    maximum bottom stress per unit area of contact (landfast ice)  
    134    REAL(wp), PUBLIC ::   rn_lfrelax       !:    relaxation time scale (s-1) to reach static friction (landfast ice)  
    135    ! 
    136    !                                     !!** ice-rheology namelist (namrhg) ** 
     131   LOGICAL , PUBLIC ::   ln_landfast_L16  !: landfast ice parameterizationfrom lemieux2016  
     132   LOGICAL , PUBLIC ::   ln_landfast_home !: landfast ice parameterizationfrom home made  
     133   REAL(wp), PUBLIC ::   rn_depfra        !:    fraction of ocean depth that ice must reach to initiate landfast ice 
     134   REAL(wp), PUBLIC ::   rn_icebfr        !:    maximum bottom stress per unit area of contact (lemieux2016) or per unit volume (home)  
     135   REAL(wp), PUBLIC ::   rn_lfrelax       !:    relaxation time scale (s-1) to reach static friction 
     136   REAL(wp), PUBLIC ::   rn_tensile       !:    isotropic tensile strength 
     137   ! 
     138   !                                     !!** ice-ridging/rafting namelist (namdyn_rdgrft) ** 
     139   REAL(wp), PUBLIC ::   rn_crhg          !: determines changes in ice strength (also used for landfast param) 
     140   ! 
     141   !                                     !!** ice-rheology namelist (namdyn_rhg) ** 
    137142   LOGICAL , PUBLIC ::   ln_aEVP          !: using adaptive EVP (T or F)  
    138143   REAL(wp), PUBLIC ::   rn_creepl        !: creep limit : has to be under 1.0e-9 
  • NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn.F90

    r10267 r10312  
    8484         WRITE(numout,*)'~~~~~~~' 
    8585      ENDIF 
    86  
    8786      !                       
    88       IF( ln_landfast ) THEN            !-- Landfast ice parameterization: define max bottom friction 
     87      IF( ln_landfast_home ) THEN      !-- Landfast ice parameterization 
    8988         tau_icebfr(:,:) = 0._wp 
    9089         DO jl = 1, jpl 
    91             WHERE( h_i(:,:,jl) > ht_n(:,:) * rn_gamma )   tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr 
     90            WHERE( h_i(:,:,jl) > ht_n(:,:) * rn_depfra )   tau_icebfr(:,:) = tau_icebfr(:,:) + a_i(:,:,jl) * rn_icebfr 
    9291         END DO 
    93          IF( iom_use('tau_icebfr') )   CALL iom_put( 'tau_icebfr', tau_icebfr )   
    94       ENDIF 
    95  
     92      ENDIF 
     93      ! 
    9694      zhmax(:,:,:) = h_i_b(:,:,:)      !-- Record max of the surrounding 9-pts ice thick. (for CALL Hbig) 
    9795      DO jl = 1, jpl 
     
    244242      !! 
    245243      NAMELIST/namdyn/ ln_dynFULL, ln_dynRHGADV, ln_dynADV, rn_uice, rn_vice,  & 
    246          &             rn_ishlat  ,                                            & 
    247          &             ln_landfast, rn_gamma , rn_icebfr, rn_lfrelax 
     244         &             rn_ishlat ,                                             & 
     245         &             ln_landfast_L16, ln_landfast_home, rn_depfra, rn_icebfr, rn_lfrelax, rn_tensile 
    248246      !!------------------------------------------------------------------- 
    249247      ! 
     
    261259         WRITE(numout,*) '~~~~~~~~~~~~' 
    262260         WRITE(numout,*) '   Namelist namdyn:' 
    263          WRITE(numout,*) '      Full ice dynamics      (rhg + adv + ridge/raft + corr)  ln_dynFULL   = ', ln_dynFULL 
    264          WRITE(numout,*) '      No ridge/raft & No cor (rhg + adv)                      ln_dynRHGADV = ', ln_dynRHGADV 
    265          WRITE(numout,*) '      Advection only         (rn_uvice + adv)                 ln_dynADV    = ', ln_dynADV 
    266          WRITE(numout,*) '         with prescribed velocity given by   (u,v)_ice = (rn_uice,rn_vice) = (', rn_uice,',', rn_vice,')' 
    267          WRITE(numout,*) '      lateral boundary condition for sea ice dynamics         rn_ishlat    = ', rn_ishlat 
    268          WRITE(numout,*) '      Landfast: param (T or F)                                ln_landfast  = ', ln_landfast 
    269          WRITE(numout,*) '         fraction of ocean depth that ice must reach          rn_gamma     = ', rn_gamma 
    270          WRITE(numout,*) '         maximum bottom stress per unit area of contact       rn_icebfr    = ', rn_icebfr 
    271          WRITE(numout,*) '         relax time scale (s-1) to reach static friction      rn_lfrelax   = ', rn_lfrelax 
     261         WRITE(numout,*) '      Full ice dynamics      (rhg + adv + ridge/raft + corr)  ln_dynFULL      = ', ln_dynFULL 
     262         WRITE(numout,*) '      No ridge/raft & No cor (rhg + adv)                      ln_dynRHGADV    = ', ln_dynRHGADV 
     263         WRITE(numout,*) '      Advection only         (rn_uvice + adv)                 ln_dynADV       = ', ln_dynADV 
     264         WRITE(numout,*) '         with prescribed velocity given by   (u,v)_ice = (rn_uice,rn_vice)    = (', rn_uice,',', rn_vice,')' 
     265         WRITE(numout,*) '      lateral boundary condition for sea ice dynamics         rn_ishlat       = ', rn_ishlat 
     266         WRITE(numout,*) '      Landfast: param from Lemieux 2016                       ln_landfast_L16 = ', ln_landfast_L16 
     267         WRITE(numout,*) '      Landfast: param from home made                          ln_landfast_home= ', ln_landfast_home 
     268         WRITE(numout,*) '         fraction of ocean depth that ice must reach          rn_depfra       = ', rn_depfra 
     269         WRITE(numout,*) '         maximum bottom stress per unit area of contact       rn_icebfr       = ', rn_icebfr 
     270         WRITE(numout,*) '         relax time scale (s-1) to reach static friction      rn_lfrelax      = ', rn_lfrelax 
     271         WRITE(numout,*) '         isotropic tensile strength                           rn_tensile      = ', rn_tensile 
    272272         WRITE(numout,*) 
    273273      ENDIF 
     
    289289      ELSEIF ( 2. < rn_ishlat                      ) THEN   ;   IF(lwp) WRITE(numout,*) '   ===>>>   ice lateral  strong-slip' 
    290290      ENDIF 
    291       !                                      !--- NO Landfast ice : set to zero once for all 
    292       IF( .NOT.ln_landfast )   tau_icebfr(:,:) = 0._wp 
     291      !                                      !--- Landfast ice 
     292      IF( .NOT.ln_landfast_L16 .AND. .NOT.ln_landfast_home )   tau_icebfr(:,:) = 0._wp 
     293      ! 
     294      IF ( ln_landfast_L16 .AND. ln_landfast_home ) THEN 
     295         CALL ctl_stop( 'ice_dyn_init: choose one and only one landfast parameterization (ln_landfast_L16 or ln_landfast_home)' ) 
     296      ENDIF 
    293297      ! 
    294298      CALL ice_dyn_rdgrft_init          ! set ice ridging/rafting parameters 
  • NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_rdgrft.F90

    r9935 r10312  
    3939 
    4040   ! Variables shared among ridging subroutines 
    41    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:) ::   closing_net     ! net rate at which area is removed    (1/s) 
    42       !                                                 ! (ridging ice area - area of new ridges) / dt 
    43    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:)   ::   opning         ! rate of opening due to divergence/shear 
    44    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:)   ::   closing_gross  ! rate at which area removed, not counting area of new ridges 
    45    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   apartf   ! participation function; fraction of ridging/closing associated w/ category n 
    46    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hrmin    ! minimum ridge thickness 
    47    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hrmax    ! maximum ridge thickness 
    48    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hraft    ! thickness of rafted ice 
    49    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   hi_hrdg  ! thickness of ridging ice / mean ridge thickness 
    50    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   aridge   ! participating ice ridging 
    51    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   araft    ! participating ice rafting 
     41   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:)     ::   closing_net     ! net rate at which area is removed    (1/s) 
     42      !                                                               ! (ridging ice area - area of new ridges) / dt 
     43   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:)     ::   opning          ! rate of opening due to divergence/shear 
     44   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:)     ::   closing_gross   ! rate at which area removed, not counting area of new ridges 
     45   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   apartf          ! participation function; fraction of ridging/closing associated w/ category n 
     46   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hrmin           ! minimum ridge thickness 
     47   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hrmax           ! maximum ridge thickness 
     48   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hraft           ! thickness of rafted ice 
     49   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hi_hrdg         ! thickness of ridging ice / mean ridge thickness 
     50   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   aridge          ! participating ice ridging 
     51   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   araft           ! participating ice rafting 
    5252   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ze_i_2d 
    5353   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   ze_s_2d 
     
    5959   LOGICAL  ::   ln_str_H79       ! ice strength parameterization (Hibler79) 
    6060   REAL(wp) ::   rn_pstar         ! determines ice strength, Hibler JPO79 
    61    REAL(wp) ::   rn_crhg          ! determines changes in ice strength 
    6261   REAL(wp) ::   rn_csrdg         ! fraction of shearing energy contributing to ridging             
    6362   LOGICAL  ::   ln_partf_lin     ! participation function linear (Thorndike et al. (1975)) 
  • NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icedyn_rhg_evp.F90

    r9935 r10312  
    118118      REAL(wp) ::   ecc2, z1_ecc2                                       ! square of yield ellipse eccenticity 
    119119      REAL(wp) ::   zalph1, z1_alph1, zalph2, z1_alph2                  ! alpha coef from Bouillon 2009 or Kimmritz 2017 
    120       REAL(wp) ::   zm1, zm2, zm3, zmassU, zmassV                       ! ice/snow mass 
     120      REAL(wp) ::   zm1, zm2, zm3, zmassU, zmassV, zvU, zvV             ! ice/snow mass and volume 
    121121      REAL(wp) ::   zdelta, zp_delf, zds2, zdt, zdt2, zdiv, zdiv2       ! temporary scalars 
    122122      REAL(wp) ::   zTauO, zTauB, zTauE, zvel                           ! temporary scalars 
     123      REAL(wp) ::   zkt                                                 ! isotropic tensile strength for landfast ice 
     124      REAL(wp) ::   zvCr                                                ! critical ice volume above which ice is landfast 
    123125      ! 
    124126      REAL(wp) ::   zresm                                               ! Maximal error on ice velocity 
     
    136138      REAL(wp), DIMENSION(jpi,jpj) ::   zmf                             ! coriolis parameter at T points 
    137139      REAL(wp), DIMENSION(jpi,jpj) ::   zTauU_ia , ztauV_ia             ! ice-atm. stress at U-V points 
     140      REAL(wp), DIMENSION(jpi,jpj) ::   zTauU_ib , ztauV_ib             ! ice-bottom stress at U-V points (landfast param) 
    138141      REAL(wp), DIMENSION(jpi,jpj) ::   zspgU , zspgV                   ! surface pressure gradient at U/V points 
    139142      REAL(wp), DIMENSION(jpi,jpj) ::   v_oceU, u_oceV, v_iceU, u_iceV  ! ocean/ice u/v component on V/U points                            
     
    255258         END DO 
    256259      END DO 
    257              
     260 
     261      ! landfast param from Lemieux(2016): add isotropic tensile strength (following Konig Beatty and Holland, 2010) 
     262      IF( ln_landfast_L16 .OR. ln_landfast_home ) THEN   ;   zkt = rn_tensile 
     263      ELSE                                               ;   zkt = 0._wp 
     264      ENDIF 
    258265      ! 
    259266      !------------------------------------------------------------------------------! 
     
    273280         zpice(:,:) = ssh_m(:,:) + ( zintn * snwice_mass(:,:) + zintb * snwice_mass_b(:,:) ) * r1_rau0 
    274281         ! 
    275       ELSE                                    !== non-embedded sea ice: use ocean surface for slope calculation ==! 
     282      ELSE                               !== non-embedded sea ice: use ocean surface for slope calculation ==! 
    276283         zpice(:,:) = ssh_m(:,:) 
    277284      ENDIF 
     
    328335      CALL lbc_lnk_multi( zmf, 'T', 1., zdt_m, 'T', 1. ) 
    329336      ! 
     337      !                                  !== Landfast ice parameterization ==! 
     338      ! 
     339      IF( ln_landfast_L16 ) THEN         !-- Lemieux 2016 
     340         DO jj = 2, jpjm1 
     341            DO ji = fs_2, fs_jpim1 
     342               ! ice thickness at U-V points 
     343               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) 
     344               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) 
     345               ! ice-bottom stress at U points 
     346               zvCr = zaU(ji,jj) * rn_depfra * hu_n(ji,jj) 
     347               zTauU_ib(ji,jj)   = rn_icebfr * MAX( 0._wp, zvU - zvCr ) * EXP( -rn_crhg * ( 1._wp - zaU(ji,jj) ) ) 
     348               ! ice-bottom stress at V points 
     349               zvCr = zaV(ji,jj) * rn_depfra * hv_n(ji,jj) 
     350               zTauV_ib(ji,jj)   = rn_icebfr * MAX( 0._wp, zvV - zvCr ) * EXP( -rn_crhg * ( 1._wp - zaV(ji,jj) ) ) 
     351               ! ice_bottom stress at T points 
     352               zvCr = at_i(ji,jj) * rn_depfra * ht_n(ji,jj) 
     353               tau_icebfr(ji,jj) = rn_icebfr * MAX( 0._wp, vt_i(ji,jj) - zvCr ) * EXP( -rn_crhg * ( 1._wp - at_i(ji,jj) ) ) 
     354            END DO 
     355         END DO 
     356         CALL lbc_lnk( tau_icebfr(:,:), 'T', 1. ) 
     357         ! 
     358      ELSEIF( ln_landfast_home ) THEN          !-- Home made 
     359         DO jj = 2, jpjm1 
     360            DO ji = fs_2, fs_jpim1 
     361               zTauU_ib(ji,jj) = tau_icebfr(ji,jj) 
     362               zTauV_ib(ji,jj) = tau_icebfr(ji,jj) 
     363            END DO 
     364         END DO 
     365         ! 
     366      ELSE                                     !-- no landfast 
     367         DO jj = 2, jpjm1 
     368            DO ji = fs_2, fs_jpim1 
     369               zTauU_ib(ji,jj) = 0._wp 
     370               zTauV_ib(ji,jj) = 0._wp 
     371            END DO 
     372         END DO 
     373      ENDIF 
     374      IF( iom_use('tau_icebfr') )   CALL iom_put( 'tau_icebfr', tau_icebfr(:,:) ) 
     375 
    330376      !------------------------------------------------------------------------------! 
    331377      ! 3) Solution of the momentum equation, iterative procedure 
     
    391437               ENDIF 
    392438                
    393                ! stress at T points 
    394                zs1(ji,jj) = ( zs1(ji,jj) * zalph1 + zp_delt(ji,jj) * ( zdiv - zdelta ) ) * z1_alph1 
    395                zs2(ji,jj) = ( zs2(ji,jj) * zalph2 + zp_delt(ji,jj) * ( zdt * z1_ecc2 ) ) * z1_alph2 
     439               ! stress at T points (zkt/=0 if landfast) 
     440               zs1(ji,jj) = ( zs1(ji,jj) * zalph1 + zp_delt(ji,jj) * ( zdiv * (1._wp + zkt) - zdelta * (1._wp - zkt) ) ) * z1_alph1 
     441               zs2(ji,jj) = ( zs2(ji,jj) * zalph2 + zp_delt(ji,jj) * ( zdt * z1_ecc2 * (1._wp + zkt) ) ) * z1_alph2 
    396442              
    397443            END DO 
     
    412458               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) ) 
    413459                
    414                ! stress at F points 
    415                zs12(ji,jj)= ( zs12(ji,jj) * zalph2 + zp_delf * ( zds(ji,jj) * z1_ecc2 ) * 0.5_wp ) * z1_alph2 
     460               ! stress at F points (zkt/=0 if landfast) 
     461               zs12(ji,jj)= ( zs12(ji,jj) * zalph2 + zp_delf * ( zds(ji,jj) * z1_ecc2 * (1._wp + zkt) ) * 0.5_wp ) * z1_alph2 
    416462 
    417463            END DO 
     
    461507                  ! 
    462508                  !                 !--- tau_bottom/v_ice 
    463                   zvel  = MAX( zepsi, SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) ) 
    464                   zTauB = - tau_icebfr(ji,jj) / zvel 
     509                  zvel  = 5.e-05_wp + SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) 
     510                  zTauB = - zTauV_ib(ji,jj) / zvel 
    465511                  ! 
    466512                  !                 !--- Coriolis at V-points (energy conserving formulation) 
     
    473519                  ! 
    474520                  !                 !--- landfast switch => 0 = static friction ; 1 = sliding friction 
    475                   rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - tau_icebfr(ji,jj) ) - SIGN( 1._wp, zTauE ) ) ) 
     521                  rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - zTauV_ib(ji,jj) ) - SIGN( 1._wp, zTauE ) ) ) 
    476522                  ! 
    477523                  IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017) 
     
    509555                  ! 
    510556                  !                 !--- tau_bottom/u_ice 
    511                   zvel  = MAX( zepsi, SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) ) 
    512                   zTauB = - tau_icebfr(ji,jj) / zvel 
     557                  zvel  = 5.e-05_wp + SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) 
     558                  zTauB = - zTauU_ib(ji,jj) / zvel 
    513559                  ! 
    514560                  !                 !--- Coriolis at U-points (energy conserving formulation) 
     
    521567                  ! 
    522568                  !                 !--- landfast switch => 0 = static friction ; 1 = sliding friction 
    523                   rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - tau_icebfr(ji,jj) ) - SIGN( 1._wp, zTauE ) ) ) 
     569                  rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - zTauU_ib(ji,jj) ) - SIGN( 1._wp, zTauE ) ) ) 
    524570                  ! 
    525571                  IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017) 
     
    559605                  ! 
    560606                  !                 !--- tau_bottom/u_ice 
    561                   zvel  = MAX( zepsi, SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) ) 
    562                   zTauB = - tau_icebfr(ji,jj) / zvel 
     607                  zvel  = 5.e-05_wp + SQRT( v_iceU(ji,jj) * v_iceU(ji,jj) + u_ice(ji,jj) * u_ice(ji,jj) ) 
     608                  zTauB = - zTauU_ib(ji,jj) / zvel 
    563609                  ! 
    564610                  !                 !--- Coriolis at U-points (energy conserving formulation) 
     
    571617                  ! 
    572618                  !                 !--- landfast switch => 0 = static friction ; 1 = sliding friction 
    573                   rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - tau_icebfr(ji,jj) ) - SIGN( 1._wp, zTauE ) ) ) 
     619                  rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, ztauE - zTauU_ib(ji,jj) ) - SIGN( 1._wp, zTauE ) ) ) 
    574620                  ! 
    575621                  IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017) 
     
    607653                  ! 
    608654                  !                 !--- tau_bottom/v_ice 
    609                   zvel  = MAX( zepsi, SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) ) 
    610                   ztauB = - tau_icebfr(ji,jj) / zvel 
     655                  zvel  = 5.e-05_wp + SQRT( v_ice(ji,jj) * v_ice(ji,jj) + u_iceV(ji,jj) * u_iceV(ji,jj) ) 
     656                  zTauB = - zTauV_ib(ji,jj) / zvel 
    611657                  ! 
    612658                  !                 !--- Coriolis at v-points (energy conserving formulation) 
     
    619665                  ! 
    620666                  !                 !--- landfast switch => 0 = static friction ; 1 = sliding friction 
    621                   rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, zTauE - tau_icebfr(ji,jj) ) - SIGN( 1._wp, zTauE ) ) ) 
     667                  rswitch = 1._wp - MIN( 1._wp, ABS( SIGN( 1._wp, zTauE - zTauV_ib(ji,jj) ) - SIGN( 1._wp, zTauE ) ) ) 
    622668                  ! 
    623669                  IF( ln_aEVP ) THEN !--- ice velocity using aEVP (Kimmritz et al 2016 & 2017) 
  • NEMO/branches/2018/dev_r9947_SI3_advection/src/ICE/icewri.F90

    r9935 r10312  
    5050      INTEGER  ::   ji, jj, jk, jl  ! dummy loop indices 
    5151      REAL(wp) ::   z2da, z2db, zrho1, zrho2 
    52       REAL(wp), DIMENSION(jpi,jpj)     ::   z2d !  2D workspace 
     52      REAL(wp), DIMENSION(jpi,jpj)     ::   z2d, zfast !  2D workspace 
    5353      REAL(wp), DIMENSION(jpi,jpj)     ::   zmsk00, zmsk05, zmsk15, zmsksn ! O%, 5% and 15% concentration mask and snow mask 
    5454      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zmsk00l, zmsksnl               ! cat masks 
     
    132132      IF( iom_use('vtau_ai'  ) )   CALL iom_put( "vtau_ai", vtau_ice * zmsk00     )   ! Wind stress term in force balance (y) 
    133133 
    134       IF( iom_use('icevel') ) THEN  
     134      IF( iom_use('icevel') .OR. iom_use('fasticepres') ) THEN  
     135        ! module of ice velocity 
    135136         DO jj = 2 , jpjm1 
    136137            DO ji = 2 , jpim1 
     
    141142         END DO 
    142143         CALL lbc_lnk( z2d, 'T', 1. ) 
    143          IF( iom_use('icevel') )   CALL iom_put( "icevel" , z2d                   )   ! ice velocity module 
     144         IF( iom_use('icevel') )   CALL iom_put( "icevel" , z2d ) 
     145 
     146        ! record presence of fast ice 
     147         WHERE( z2d(:,:) < 5.e-04_wp .AND. zmsk00(:,:) == 1._wp ) ; zfast(:,:) = 1._wp 
     148         ELSEWHERE                                                ; zfast(:,:) = 0._wp 
     149         END WHERE 
     150         IF( iom_use('fasticepres') )   CALL iom_put( "fasticepres" , zfast ) 
    144151      ENDIF 
    145152 
Note: See TracChangeset for help on using the changeset viewer.