Changeset 12261


Ignore:
Timestamp:
2019-12-17T14:43:09+01:00 (10 months ago)
Author:
stefryn
Message:

add yield surface output and include effect on ridging

Location:
NEMO/branches/2019/dev_r11842_SI3-10_EAP
Files:
8 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2019/dev_r11842_SI3-10_EAP/cfgs/SHARED/field_def_nemo-ice.xml

    r11878 r12261  
    8282          <field id="iceshe"       long_name="Maximum shear of sea-ice velocity field"                 standard_name="maximum_shear_of_sea_ice_velocity"         unit="s-1"  /> 
    8383          <field id="aniso"        long_name="anisotropy of sea ice floe orientation (0.5 - 1)"        standard_name="anisotropy"                                unit=""     /> 
     84          <field id="yield11"      long_name="yield surface tensor component 11"                       standard_name="yield11"                                   unit="N/m"  /> 
     85          <field id="yield22"      long_name="yield surface tensor component 22"                       standard_name="yield22"                                   unit="N/m"  /> 
     86          <field id="yield12"      long_name="yield surface tensor component 12"                       standard_name="yield12"                                   unit="N/m"  /> 
    8487      
    8588     <!-- surface heat fluxes --> 
     
    416419          <field field_ref="vfxthin"          name="vfxthin" /> 
    417420    
    418      <!-- diag error for negative ice volume after advection --> 
    419      <field field_ref="iceneg_pres"      name="sineg_pres" /> 
    420      <field field_ref="iceneg_volu"      name="sineg_volu" /> 
    421      <field field_ref="iceneg_hfx"       name="sineg_hfx"  /> 
    422421        </field_group> 
    423422 
  • NEMO/branches/2019/dev_r11842_SI3-10_EAP/src/ICE/ice.F90

    r11951 r12261  
    146146   ! 
    147147   !                                     !!** ice-rheology namelist (namdyn_rhg) ** 
     148   LOGICAL , PUBLIC ::   ln_rhg_EVP       ! EVP rheology switch, used for rdgrft and rheology 
     149   LOGICAL , PUBLIC ::   ln_rhg_EAP       ! EAP rheology switch, used for rdgrft and rheology 
    148150   LOGICAL , PUBLIC ::   ln_aEVP          !: using adaptive EVP (T or F)  
    149151   REAL(wp), PUBLIC ::   rn_creepl        !: creep limit : has to be under 1.0e-9 
  • NEMO/branches/2019/dev_r11842_SI3-10_EAP/src/ICE/icedyn_rdgrft.F90

    r11732 r12261  
    2323   USE icevar         ! sea-ice: operations 
    2424   USE icectl         ! sea-ice: control prints 
     25   !USE icedyn_rhg     ! sea-ice: rheology choice 
    2526   ! 
    2627   USE in_out_manager ! I/O manager 
     
    138139      INTEGER , DIMENSION(jpij) ::   iptidx        ! compute ridge/raft or not 
    139140      REAL(wp), DIMENSION(jpij) ::   zdivu, zdelt  ! 1D divu_i & delta_i 
     141      REAL(wp), DIMENSION(jpij) ::   conv          ! 1D rdg_conv (if EAP rheology) 
    140142      ! 
    141143      INTEGER, PARAMETER ::   jp_itermax = 20     
     
    174176         
    175177         ! just needed here 
    176          CALL tab_2d_1d( npti, nptidx(1:npti), zdelt   (1:npti)      , delta_i ) 
     178         CALL tab_2d_1d( npti, nptidx(1:npti), zdelt   (1:npti)      , delta_i  ) 
     179         CALL tab_2d_1d( npti, nptidx(1:npti), conv    (1:npti)      , rdg_conv ) 
    177180         ! needed here and in the iteration loop 
    178181         CALL tab_2d_1d( npti, nptidx(1:npti), zdivu   (1:npti)      , divu_i) ! zdivu is used as a work array here (no change in divu_i) 
     
    185188            !                                                        - ice area added in new ridges 
    186189            closing_net(ji) = rn_csrdg * 0.5_wp * ( zdelt(ji) - ABS( zdivu(ji) ) ) - MIN( zdivu(ji), 0._wp ) 
     190            IF( ln_rhg_EVP )  closing_net(ji) = rn_csrdg * 0.5_wp * ( zdelt(ji) - ABS( zdivu(ji) ) ) - MIN( zdivu(ji), 0._wp ) 
     191            IF( ln_rhg_EAP )  closing_net(ji) = conv(ji) 
    187192            ! 
    188193            IF( zdivu(ji) < 0._wp )   closing_net(ji) = MAX( closing_net(ji), -zdivu(ji) )   ! make sure the closing rate is large enough 
  • NEMO/branches/2019/dev_r11842_SI3-10_EAP/src/ICE/icedyn_rhg.F90

    r11877 r12261  
    3737 
    3838   ! ** namelist (namrhg) ** 
    39    LOGICAL ::   ln_rhg_EVP       ! EVP rheology 
    40    LOGICAL ::   ln_rhg_EAP       ! EVP rheology 
     39!   LOGICAL ::   ln_rhg_EVP       ! EVP rheology, used by icedyn_rdgrft 
     40!   LOGICAL ::   ln_rhg_EAP       ! EAP rheology, used by icedyn_rdgrft 
     41   !LOGICAL, PUBLIC ::   ln_rhg_EVP       ! EVP rheology, used by icedyn_rdgrft 
     42   !LOGICAL, PUBLIC ::   ln_rhg_EAP       ! EAP rheology, used by icedyn_rdgrft 
    4143   ! 
    4244   !! * Substitutions 
  • NEMO/branches/2019/dev_r11842_SI3-10_EAP/src/ICE/icedyn_rhg_eap.F90

    r11951 r12261  
    145145      ! 
    146146      REAL(wp), DIMENSION(jpi,jpj) ::   stress12tmp                     ! anisotropic stress tensor component for regridding 
     147      REAL(wp), DIMENSION(jpi,jpj) ::   yield11, yield22, yield12       ! yield surface tensor for history 
    147148      REAL(wp), DIMENSION(jpi,jpj) ::   zp_delt                         ! P/delta at T points 
    148149      REAL(wp), DIMENSION(jpi,jpj) ::   zbeta                           ! beta coef from Kimmritz 2017 
     
    384385                  &   + zds(ji,jj-1) * e1e2f(ji,jj-1) + zds(ji-1,jj-1) * e1e2f(ji-1,jj-1)  & 
    385386                  &   ) * 0.25_wp * r1_e1e2t(ji,jj) 
     387               !WRITE(numout,*) 'shear output', ji, jj, zdsT 
     388                
    386389 
    387390               ! shear**2 at T points (doc eq. A16) 
     
    409412                                       alphar, alphas) 
    410413               IF (jter == nn_nevp) THEN 
     414                  yield11(ji,jj) = 0.5_wp * (stressptmp + stressmtmp) 
     415                  yield22(ji,jj) = 0.5_wp * (stressptmp - stressmtmp) 
     416                  yield12(ji,jj) = stress12tmp(ji,jj) 
    411417                  rdg_conv(ji,jj) = -min( alphar, 0._wp)     
    412418               ENDIF 
     
    421427               !   gamma = 0.5*P/(delta+creepl) * (c*pi)**2/Area * dt/m 
    422428               !   alpha = beta = sqrt(4*gamma) 
    423                !IF( ln_aEVP ) THEN 
    424                !   zalph1   = MAX( 50._wp, rpi * SQRT( 0.5_wp * zp_delt(ji,jj) * r1_e1e2t(ji,jj) * zdt_m(ji,jj) ) ) 
    425                !   z1_alph1 = 1._wp / ( zalph1 + 1._wp ) 
    426                !ENDIF 
     429               IF( ln_aEVP ) THEN 
     430                  zalph1   = MAX( 50._wp, rpi * SQRT( 0.5_wp * zp_delt(ji,jj) * r1_e1e2t(ji,jj) * zdt_m(ji,jj) ) ) 
     431                  z1_alph1 = 1._wp / ( zalph1 + 1._wp ) 
     432               ENDIF 
    427433                
    428434               ! stress at T points (zkt/=0 if landfast) 
     
    437443         END DO 
    438444         !CALL lbc_lnk( 'icedyn_rhg_eap', zp_delt, 'T', 1. ) 
    439          CALL lbc_lnk( 'icedyn_rhg_eap', stress12tmp, 'T', 1. ) 
     445         CALL lbc_lnk_multi( 'icedyn_rhg_eap', stress12tmp, 'T', 1., rdg_conv, 'T', 1. ) 
    440446 
    441447         DO jj = 1, jpjm1 
     
    447453 
    448454               ! alpha & beta for aEVP 
    449                !IF( ln_aEVP ) THEN 
    450                !   zalph1   = MAX( 50._wp, rpi * SQRT( 0.5_wp * zp_delt(ji,jj) * r1_e1e2t(ji,jj) * zdt_m(ji,jj) ) ) 
    451                !   z1_alph1 = 1._wp / ( zalph1 + 1._wp ) 
    452                !ENDIF 
     455               IF( ln_aEVP ) THEN 
     456                  zalph1   = MAX( 50._wp, rpi * SQRT( 0.5_wp * zp_delt(ji,jj) * r1_e1e2t(ji,jj) * zdt_m(ji,jj) ) ) 
     457                  z1_alph1 = 1._wp / ( zalph1 + 1._wp ) 
     458               ENDIF 
    453459                
    454460               ! stress at F points (zkt/=0 if landfast) 
     
    828834 
    829835         DEALLOCATE( zsig1 , zsig2 , zsig3 ) 
     836      ENDIF 
     837 
     838      ! --- yieldcurve --- ! 
     839      IF( iom_use('yield11') .OR. iom_use('yield12') .OR. iom_use('yield22')) THEN 
     840         CALL iom_put( 'yield11', yield11 * zmsk00 ) 
     841         CALL iom_put( 'yield22', yield22 * zmsk00 ) 
     842         CALL iom_put( 'yield12', yield12 * zmsk00 ) 
    830843      ENDIF 
    831844 
  • NEMO/branches/2019/dev_r11842_SI3-10_EAP/tests/ICE_ADV2D/EXPREF/file_def_nemo-ice.xml

    r10516 r12261  
    2323        <field field_ref="icethic"          name="sithic_max" operation="maximum" /> 
    2424        <field field_ref="icethic"          name="sithic_min" operation="minimum" /> 
    25         <field field_ref="iceneg_pres"      name="sineg_pres" /> 
    26         <field field_ref="iceneg_volu"      name="sineg_volu" /> 
    2725        <field field_ref="fasticepres"      name="fasticepres" /> 
    2826        <field field_ref="icevolu"          name="sivolu" /> 
  • NEMO/branches/2019/dev_r11842_SI3-10_EAP/tests/ICE_ADV2D/EXPREF/namelist_ice_cfg

    r10535 r12261  
    5050&namdyn_adv     !   Ice advection 
    5151!------------------------------------------------------------------------------ 
    52    ln_adv_Pra       = .false.         !  Advection scheme (Prather) 
    53    ln_adv_UMx       = .true.          !  Advection scheme (Ultimate-Macho) 
     52   ln_adv_Pra       = .true.         !  Advection scheme (Prather) 
     53   ln_adv_UMx       = .false.          !  Advection scheme (Ultimate-Macho) 
    5454      nn_UMx        =   5             !     order of the scheme for UMx (1-5 ; 20=centered 2nd order) 
    5555/ 
  • NEMO/branches/2019/dev_r11842_SI3-10_EAP/tests/ICE_ADV2D/MY_SRC/usrdef_nam.F90

    r11536 r12261  
    1414   !!   usr_def_hgr   : initialize the horizontal mesh  
    1515   !!---------------------------------------------------------------------- 
    16    USE dom_oce  , ONLY: nimpp , njmpp            ! i- & j-indices of the local domain 
     16   USE dom_oce  , ONLY: nimpp , njmpp , Agrif_Root           ! i- & j-indices of the local domain 
    1717   USE par_oce        ! ocean space and time domain 
    1818   USE phycst         ! physical constants 
Note: See TracChangeset for help on using the changeset viewer.