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 15542 for NEMO/branches/UKMO – NEMO

Changeset 15542 for NEMO/branches/UKMO


Ignore:
Timestamp:
2021-11-26T15:51:22+01:00 (2 years ago)
Author:
emmafiedler
Message:

Exponential ridge ice redistribution

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/UKMO/NEMO_4.0.4_ice_strength/src/ICE/icedyn_rdgrft.F90

    r15507 r15542  
    4949   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hrmin           ! minimum ridge thickness 
    5050   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hrmax           ! maximum ridge thickness 
     51   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hrexp           ! e-folding ridge thickness 
    5152   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hraft           ! thickness of rafted ice 
    5253   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   hi_hrdg         ! thickness of ridging ice / mean ridge thickness 
     
    6465   LOGICAL  ::   ln_str_R75       ! ice strength parameterization (Rothrock75) 
    6566   REAL(wp) ::   rn_Cf            !    ratio of ridging work to PE change in ridging, R75 
    66    REAL(wp) ::   rn_murdg         !    gives e-folding scale of ridged ice (m^.5), R75 
     67   INTEGER  ::   rdg_redistr      ! redistribution of ridged ice (0 = uniform (Hibler 1980), 1 = exponential) 
     68   REAL(wp) ::   rn_murdg         !    gives e-folding scale of ridged ice (m^.5) for rdg_redistr = 1 
    6769   REAL(wp) ::   rn_csrdg         ! fraction of shearing energy contributing to ridging             
    6870   LOGICAL  ::   ln_partf_lin     ! participation function linear (Thorndike et al. (1975)) 
     
    9496      ALLOCATE( closing_net(jpij)  , opning(jpij)      , closing_gross(jpij) ,  & 
    9597         &      apartf(jpij,0:jpl) , hrmin(jpij,jpl)   , hrmax(jpij,jpl)     , zhi(jpij,jpl),  & 
    96          &      zaksum(jpij)       , hraft(jpij,jpl)   , & 
     98         &      zaksum(jpij)       , hraft(jpij,jpl)   , hrexp(jpij,jpl)     , & 
    9799         &      hi_hrdg(jpij,jpl)  , araft(jpij,jpl)   , aridge(jpij,jpl)    , zasum(jpij), & 
    98100         &      ze_i_2d(jpij,nlay_i,jpl), ze_s_2d(jpij,nlay_s,jpl), STAT=ice_dyn_rdgrft_alloc ) 
     
    450452               zhmean         = MAX( SQRT( rn_hstar * zhi(ji,jl) ), zhi(ji,jl) * hrdg_hi_min ) 
    451453               hrmin  (ji,jl) = MIN( 2._wp * zhi(ji,jl), 0.5_wp * ( zhmean + zhi(ji,jl) ) ) 
    452                hrmax  (ji,jl) = 2._wp * zhmean - hrmin(ji,jl) 
    453454               hraft  (ji,jl) = zhi(ji,jl) * zfac 
    454                hi_hrdg(ji,jl) = zhi(ji,jl) / zhmean 
     455               ! 
     456               IF (rdg_redistr == 0) THEN 
     457                 hrmax  (ji,jl) = 2._wp * zhmean - hrmin(ji,jl) 
     458                 hi_hrdg(ji,jl) = zhi(ji,jl) / zhmean 
     459               ELSE IF (rdg_redistr == 1) THEN 
     460                 hrexp  (ji,jl) = rn_murdg * SQRT( zhi(ji,jl) )                
     461                 hi_hrdg(ji,jl) = zhi(ji,jl) / ( hrmin(ji,jl) + hrexp(ji,jl) )      
     462               END IF     
    455463               ! 
    456464               ! Normalization factor : zaksum, ensures mass conservation 
     
    461469               hrmax  (ji,jl) = 0._wp  
    462470               hraft  (ji,jl) = 0._wp  
     471               hrexp  (ji,jl) = 0._wp 
    463472               hi_hrdg(ji,jl) = 1._wp 
    464473               zaksum (ji)    = 0._wp 
     
    516525      INTEGER  ::   ji, jj, jl, jl1, jl2, jk   ! dummy loop indices 
    517526      REAL(wp) ::   hL, hR, farea              ! left and right limits of integration and new area going to jl2 
     527      REAL(wp) ::   expL, expR                 ! exponentials involving hL, hR 
    518528      REAL(wp) ::   vsw                        ! vol of water trapped into ridges 
    519529      REAL(wp) ::   afrdg, afrft               ! fraction of category area ridged/rafted  
     
    684694         itest_rdg(1:npti) = 0 
    685695         itest_rft(1:npti) = 0 
     696         ! 
    686697         DO jl2  = 1, jpl  
    687698            ! 
    688699            DO ji = 1, npti 
    689700 
    690                IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN 
    691  
    692                   ! Compute the fraction of ridged ice area and volume going to thickness category jl2 
     701              IF( apartf(ji,jl1) > 0._wp .AND. closing_gross(ji) > 0._wp ) THEN 
     702                 ! 
     703                 ! Compute the fraction of ridged ice area and volume going to thickness category jl2 
     704                 ! 
     705                 IF (rdg_redistr == 0) THEN ! Hibler 1980 uniform formulation 
     706                 ! 
    693707                  IF( hrmin(ji,jl1) <= hi_max(jl2) .AND. hrmax(ji,jl1) > hi_max(jl2-1) ) THEN 
    694708                     hL = MAX( hrmin(ji,jl1), hi_max(jl2-1) ) 
     
    702716                     fvol(ji) = 0._wp                   
    703717                  ENDIF 
     718                  ! 
     719                 ELSE IF (rdg_redistr == 1) THEN ! exponential formulation 
     720                  ! 
     721                  IF (jl2 < jpl) THEN 
     722                  ! 
     723                    IF (hrmin(ji,jl1) <= hi_max(jl2)) THEN 
     724                       hL       = MAX( hrmin(ji,jl1), hi_max(jl2-1) ) 
     725                       hR       = hi_max(jl2) 
     726                       expL     = EXP( -( hL - hrmin(ji,jl1) ) / hrexp(ji,jl1) ) 
     727                       expR     = EXP( -( hR - hrmin(ji,jl1) ) / hrexp(ji,jl1) ) 
     728                       farea    = expL - expR 
     729                       fvol(ji) = ( ( hL + hrexp(ji,jl1) ) * expL  & 
     730                                  - ( hR + hrexp(ji,jl1) ) * expR ) / ( hrmin(ji,jl1) + hrexp(ji,jl1) ) 
     731                     ELSE 
     732                       farea    = 0._wp  
     733                       fvol(ji) = 0._wp   
     734                     END IF              
     735                     !                  
     736                  ELSE             ! jl2 = jpl 
     737                    ! 
     738                    hL       = MAX( hrmin(ji,jl1), hi_max(jl2-1) ) 
     739                    expL     = EXP(-( hL - hrmin(ji,jl1) ) / hrexp(ji,jl1) ) 
     740                    farea    = expL 
     741                    fvol(ji) = ( hL + hrexp(ji,jl1) ) * expL / ( hrmin(ji,jl1) + hrexp(ji,jl1) )   
     742                    ! 
     743                  END IF            ! jl2 < jpl 
     744                  !  
     745                 END IF             ! rdg_redistr 
    704746 
    705747                  ! Compute the fraction of rafted ice area and volume going to thickness category jl2 
     
    848890            DO ji = 1, npti           ! grid cells 
    849891              IF ( apartf(ji,jl) > 0._wp ) THEN 
    850  
    851                 ! Uniform distribution of ridged ice 
    852                 h2rdg = (1._wp/3._wp) * (hrmax(ji,jl)**3._wp - hrmin(ji,jl)**3._wp)  & 
    853                          / (hrmax(ji,jl) - hrmin(ji,jl)) 
    854  
     892                ! 
     893                IF (rdg_redistr == 0) THEN   ! Uniform distribution of ridged ice                
     894                  h2rdg = (1._wp/3._wp) * (hrmax(ji,jl)**3._wp - hrmin(ji,jl)**3._wp)  & 
     895                           / (hrmax(ji,jl) - hrmin(ji,jl)) 
     896                ! 
     897                ELSE IF (rdg_redistr == 1) THEN   ! Exponential distribution of ridged ice 
     898        h2rdg = hrmin(ji,jl) * hrmin(ji,jl) &  
     899                + 2._wp * hrmin(ji,jl) * hrexp(ji,jl) &  
     900                + 2._wp * hrexp(ji,jl) * hrexp(ji,jl)  
     901                END IF 
     902                 
    855903                dh2rdg = -zhi(ji,jl)*zhi(ji,jl) + h2rdg * hi_hrdg(ji,jl) 
    856904    
     
    10391087      !! 
    10401088      NAMELIST/namdyn_rdgrft/ ln_str_H79, rn_pstar,          & 
    1041          &                    ln_str_R75, rn_Cf, rn_murdg,   & 
     1089         &                    ln_str_R75, rn_Cf,             & 
     1090         &                    rdg_redistr, rn_murdg,         & 
    10421091         &                    rn_crhg, rn_csrdg,             & 
    10431092         &                    ln_partf_lin, rn_gstar,        & 
     
    10651114         WRITE(numout,*) '      ice strength parameterization Rothrock (1975)            ln_str_R75   = ', ln_str_R75 
    10661115         WRITE(numout,*) '            ratio of ridging work to PE change in ridging      rn_Cf        = ', rn_Cf 
    1067          WRITE(numout,*) '            gives e-folding scale of ridged ice (m^.5)         rn_murdg     = ', rn_murdg 
     1116         WRITE(numout,*) '      redistribution of ridged ice                             rdg_redistr  = ', rdg_redistr 
     1117         WRITE(numout,*) '            (0 = uniform (Hibler 1980), 1 = exponential)                      ' 
     1118         WRITE(numout,*) '            e-folding scale of ridged ice for rdg_redistr=1    rn_murdg     = ', rn_murdg 
    10681119         WRITE(numout,*) '      Fraction of shear energy contributing to ridging         rn_csrdg     = ', rn_csrdg  
    10691120         WRITE(numout,*) '      linear ridging participation function                    ln_partf_lin = ', ln_partf_lin 
     
    10911142      ENDIF 
    10921143      ! 
     1144!      IF ( ( rdg_redistr .NE. 0 ) .OR. ( rdg_redistr .NE. 1 ) ) THEN  ! why doesn't this work? 
     1145!         CALL ctl_stop( 'ice_dyn_rdgrft_init: ice redistribution function rdg_redistr must be 0 (uniform) or 1 (exponential)' ) 
     1146!      ENDIF 
     1147      ! 
    10931148      IF( .NOT. ln_icethd ) THEN 
    10941149         rn_porordg = 0._wp 
Note: See TracChangeset for help on using the changeset viewer.