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.
Changeset 3791 – NEMO

Changeset 3791


Ignore:
Timestamp:
2013-02-10T14:16:41+01:00 (11 years ago)
Author:
gm
Message:

dev_MERGE_2012: #1041 : corrected diagnostic for the ice shear rate

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2012/dev_MERGE_2012/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90

    r3680 r3791  
    123123      REAL(wp) ::   z0, zr, zcca, zccb ! temporary scalars 
    124124      REAL(wp) ::   zu_ice2, zv_ice1   ! 
    125       REAL(wp) ::   zddc, zdtc, zdst   ! delta on corners and on centre 
     125      REAL(wp) ::   zddc, zdtc, zzdst   ! delta on corners and on centre 
    126126      REAL(wp) ::   zdsshx, zdsshy     ! term for the gradient of ocean surface 
    127127      REAL(wp) ::   sigma1, sigma2     ! internal ice stress 
     
    147147      REAL(wp), POINTER, DIMENSION(:,:) ::   zdd   , zdt      ! Divergence and tension at centre of grid cells 
    148148      REAL(wp), POINTER, DIMENSION(:,:) ::   zds              ! Shear on northeast corner of grid cells 
     149      REAL(wp), POINTER, DIMENSION(:,:) ::   zdst             ! Shear on centre of grid cells 
    149150      REAL(wp), POINTER, DIMENSION(:,:) ::   deltat, deltac   ! Delta at centre and corners of grid cells 
    150151      REAL(wp), POINTER, DIMENSION(:,:) ::   zs1   , zs2      ! Diagonal stress tensor components zs1 and zs2  
     
    153154      REAL(wp), POINTER, DIMENSION(:,:) ::   zpice            ! array used for the calculation of ice surface slope: 
    154155                                                              !   ocean surface (ssh_m) if ice is not embedded 
    155                                                               !   ice top surface if ice is embedded 
    156        
     156                                                              !   ice top surface if ice is embedded    
    157157      !!------------------------------------------------------------------- 
    158158 
    159159      CALL wrk_alloc( jpi,jpj, zpresh, zfrld1, zmass1, zcorl1, za1ct , zpreshc, zfrld2, zmass2, zcorl2, za2ct ) 
    160160      CALL wrk_alloc( jpi,jpj, zc1   , u_oce1, u_oce2, u_ice2, zusw  , v_oce1 , v_oce2, v_ice1                ) 
    161       CALL wrk_alloc( jpi,jpj, zf1   , deltat, zu_ice, zf2   , deltac, zv_ice , zdd   , zdt    , zds          ) 
     161      CALL wrk_alloc( jpi,jpj, zf1   , deltat, zu_ice, zf2   , deltac, zv_ice , zdd   , zdt    , zds  , zdst  ) 
    162162      CALL wrk_alloc( jpi,jpj, zdd   , zdt   , zds   , zs1   , zs2   , zs12   , zresr , zpice                 ) 
    163163 
     
    405405 
    406406               !- Calculate Delta at centre of grid cells 
    407                zdst       = (  e2u(ji  , jj) * v_ice1(ji  ,jj)          & 
     407               zzdst      = (  e2u(ji  , jj) * v_ice1(ji  ,jj)          & 
    408408                  &          - e2u(ji-1, jj) * v_ice1(ji-1,jj)          & 
    409409                  &          + e1v(ji, jj  ) * u_ice2(ji,jj  )          & 
     
    412412                  &         / area(ji,jj) 
    413413 
    414                delta = SQRT( zdd(ji,jj)*zdd(ji,jj) + ( zdt(ji,jj)*zdt(ji,jj) + zdst*zdst ) * usecc2 )   
    415                deltat(ji,jj) = MAX( SQRT(zdd(ji,jj)**2 + (zdt(ji,jj)**2 + zdst**2)*usecc2), creepl ) 
     414               delta = SQRT( zdd(ji,jj)*zdd(ji,jj) + ( zdt(ji,jj)*zdt(ji,jj) + zzdst*zzdst ) * usecc2 )   
     415               deltat(ji,jj) = MAX( SQRT(zdd(ji,jj)**2 + (zdt(ji,jj)**2 + zzdst**2)*usecc2), creepl ) 
     416!!gm faster to replace the line above with simply: 
     417!!                deltat(ji,jj) = MAX( delta, creepl ) 
     418!!gm end  
    416419 
    417420               !-Calculate stress tensor components zs1 and zs2  
     
    711714                  &        * tmi(ji+1,jj) * tmi(ji+1,jj+1) 
    712715 
    713                zdst       = (  e2u( ji  , jj   ) * v_ice1(ji,jj)          & 
    714                   &          - e2u( ji-1, jj   ) * v_ice1(ji-1,jj)         & 
    715                   &          + e1v( ji  , jj   ) * u_ice2(ji,jj)           & 
    716                   &          - e1v( ji  , jj-1 ) * u_ice2(ji,jj-1)         &  
    717                   &          )                                             & 
    718                   &         / area(ji,jj) 
    719  
    720                deltat(ji,jj) = SQRT( zdd(ji,jj)*zdd(ji,jj) + &  
    721                   &                          ( zdt(ji,jj)*zdt(ji,jj) + zdst*zdst ) * usecc2 &  
     716               zdst(ji,jj) = (  e2u( ji  , jj   ) * v_ice1(ji  ,jj  )    & 
     717                  &           - e2u( ji-1, jj   ) * v_ice1(ji-1,jj  )    & 
     718                  &           + e1v( ji  , jj   ) * u_ice2(ji  ,jj  )    & 
     719                  &           - e1v( ji  , jj-1 ) * u_ice2(ji  ,jj-1)  ) / area(ji,jj) 
     720 
     721               deltat(ji,jj) = SQRT(    zdd(ji,jj)*zdd(ji,jj)   &  
     722                  &                 + ( zdt(ji,jj)*zdt(ji,jj) + zdst(ji,jj)*zdst(ji,jj) ) * usecc2 &  
    722723                  &                          ) + creepl 
    723724 
     
    732733      ! 
    733734      ! * Invariants of the stress tensor are required for limitd_me 
    734       ! accelerates convergence and improves stability 
     735      !   (accelerates convergence and improves stability) 
    735736      DO jj = k_j1+1, k_jpj-1 
    736737         DO ji = fs_2, fs_jpim1 
    737738            divu_i (ji,jj) = zdd   (ji,jj) 
    738739            delta_i(ji,jj) = deltat(ji,jj) 
    739             shear_i(ji,jj) = zds   (ji,jj) 
     740            shear_i(ji,jj) = SQRT( zdt(ji,jj) * zdt(ji,jj) + zdst(ji,jj) * zdst(ji,jj) ) 
    740741         END DO 
    741742      END DO 
    742  
    743       ! Lateral boundary condition 
    744       CALL lbc_lnk( divu_i (:,:), 'T', 1. ) 
     743      CALL lbc_lnk( divu_i (:,:), 'T', 1. )      ! Lateral boundary condition 
    745744      CALL lbc_lnk( delta_i(:,:), 'T', 1. ) 
    746       CALL lbc_lnk( shear_i(:,:), 'F', 1. ) 
     745      CALL lbc_lnk( shear_i(:,:), 'T', 1. ) 
    747746 
    748747      ! * Store the stress tensor for the next time step 
     
    790789      CALL wrk_dealloc( jpi,jpj, zpresh, zfrld1, zmass1, zcorl1, za1ct , zpreshc, zfrld2, zmass2, zcorl2, za2ct ) 
    791790      CALL wrk_dealloc( jpi,jpj, zc1   , u_oce1, u_oce2, u_ice2, zusw  , v_oce1 , v_oce2, v_ice1                ) 
    792       CALL wrk_dealloc( jpi,jpj, zf1   , deltat, zu_ice, zf2   , deltac, zv_ice , zdd   , zdt    , zds          ) 
     791      CALL wrk_dealloc( jpi,jpj, zf1   , deltat, zu_ice, zf2   , deltac, zv_ice , zdd   , zdt    , zds  , zdst  ) 
    793792      CALL wrk_dealloc( jpi,jpj, zdd   , zdt   , zds   , zs1   , zs2   , zs12   , zresr , zpice                 ) 
    794793 
Note: See TracChangeset for help on using the changeset viewer.