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 5080 for branches/2015/dev_r5044_CNRS_LIM3CLEAN/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90 – NEMO

Ignore:
Timestamp:
2015-02-11T18:27:52+01:00 (9 years ago)
Author:
clem
Message:

LIM3 cleaning again (almost done)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2015/dev_r5044_CNRS_LIM3CLEAN/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90

    r5078 r5080  
    123123      !!  and Elizabeth C. Hunke, LANL are gratefully acknowledged 
    124124      !!--------------------------------------------------------------------! 
    125       INTEGER ::   ji, jj, jk, jl   ! dummy loop index 
    126       INTEGER ::   niter, nitermax = 20   ! local integer  
    127       LOGICAL  ::   asum_error            ! flag for asum .ne. 1 
     125      INTEGER  ::   ji, jj, jk, jl        ! dummy loop index 
     126      INTEGER  ::   niter                 ! local integer  
    128127      INTEGER  ::   iterate_ridging       ! if true, repeat the ridging 
    129       REAL(wp) ::   w1, tmpfac            ! local scalar 
     128      REAL(wp) ::   za, zfac              ! local scalar 
    130129      CHARACTER (len = 15) ::   fieldid 
    131130      REAL(wp), POINTER, DIMENSION(:,:) ::   closing_net     ! net rate at which area is removed    (1/s) 
     
    137136      REAL(wp), POINTER, DIMENSION(:,:) ::   esnow_mlt       ! energy needed to melt snow in ocean (J m-2) 
    138137      REAL(wp), POINTER, DIMENSION(:,:) ::   vt_i_init, vt_i_final  !  ice volume summed over categories 
     138      ! 
     139      INTEGER, PARAMETER ::   nitermax = 20     
    139140      ! 
    140141      REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b  
     
    236237               ! would be removed.  Reduce the opening rate proportionately. 
    237238               IF ( ato_i(ji,jj) > epsi10 .AND. athorn(ji,jj,0) > 0.0 ) THEN 
    238                   w1 = athorn(ji,jj,0) * closing_gross(ji,jj) * rdt_ice 
    239                   IF ( w1 > ato_i(ji,jj)) THEN 
    240                      tmpfac = ato_i(ji,jj) / w1 
    241                      closing_gross(ji,jj) = closing_gross(ji,jj) * tmpfac 
    242                      opning(ji,jj) = opning(ji,jj) * tmpfac 
    243                   ENDIF !w1 
    244                ENDIF !at0i and athorn 
    245  
    246             END DO ! ji 
    247          END DO ! jj 
     239                  za = athorn(ji,jj,0) * closing_gross(ji,jj) * rdt_ice 
     240                  IF ( za > ato_i(ji,jj)) THEN 
     241                     zfac = ato_i(ji,jj) / za 
     242                     closing_gross(ji,jj) = closing_gross(ji,jj) * zfac 
     243                     opning(ji,jj) = opning(ji,jj) * zfac 
     244                  ENDIF 
     245               ENDIF 
     246 
     247            END DO 
     248         END DO 
    248249 
    249250         ! correction to closing rate / opening if excessive ice removal 
     
    256257               DO ji = 1, jpi 
    257258                  IF ( a_i(ji,jj,jl) > epsi10 .AND. athorn(ji,jj,jl) > 0._wp )THEN 
    258                      w1 = athorn(ji,jj,jl) * closing_gross(ji,jj) * rdt_ice 
    259                      IF ( w1  >  a_i(ji,jj,jl) ) THEN 
    260                         tmpfac = a_i(ji,jj,jl) / w1 
    261                         closing_gross(ji,jj) = closing_gross(ji,jj) * tmpfac 
    262                         opning       (ji,jj) = opning       (ji,jj) * tmpfac 
     259                     za = athorn(ji,jj,jl) * closing_gross(ji,jj) * rdt_ice 
     260                     IF ( za  >  a_i(ji,jj,jl) ) THEN 
     261                        zfac = a_i(ji,jj,jl) / za 
     262                        closing_gross(ji,jj) = closing_gross(ji,jj) * zfac 
     263                        opning       (ji,jj) = opning       (ji,jj) * zfac 
    263264                     ENDIF 
    264265                  ENDIF 
    265                END DO !ji 
    266             END DO ! jj 
    267          END DO !jl 
     266               END DO 
     267            END DO 
     268         END DO 
    268269 
    269270         ! 3.3 Redistribute area, volume, and energy. 
     
    322323      ! Convert ridging rate diagnostics to correct units. 
    323324      ! Update fresh water and heat fluxes due to snow melt. 
    324  
    325       asum_error = .false.  
    326  
    327325      DO jj = 1, jpj 
    328326         DO ji = 1, jpi 
    329  
    330             IF(ABS(asum(ji,jj) - kamax) > epsi10 ) asum_error = .true. 
    331327 
    332328            dardg1dt(ji,jj) = dardg1dt(ji,jj) * r1_rdtice 
     
    434430      !!---------------------------------------------------------------------- 
    435431      INTEGER, INTENT(in) ::   kstrngth    ! = 1 for Rothrock formulation, 0 for Hibler (1979) 
    436  
    437       INTEGER ::   ji,jj, jl   ! dummy loop indices 
    438       INTEGER ::   ksmooth     ! smoothing the resistance to deformation 
    439       INTEGER ::   numts_rm    ! number of time steps for the P smoothing 
    440       REAL(wp) ::   hi, zw1, zp, zdummy, zzc, z1_3   ! local scalars 
     432      INTEGER             ::   ji,jj, jl   ! dummy loop indices 
     433      INTEGER             ::   ksmooth     ! smoothing the resistance to deformation 
     434      INTEGER             ::   numts_rm    ! number of time steps for the P smoothing 
     435      REAL(wp)            ::   zhi, zp, z1_3  ! local scalars 
    441436      REAL(wp), POINTER, DIMENSION(:,:) ::   zworka   ! temporary array used here 
    442437      !!---------------------------------------------------------------------- 
     
    464459                  ! 
    465460                  IF( a_i(ji,jj,jl) > epsi10 .AND. athorn(ji,jj,jl) > 0._wp ) THEN 
    466                      hi = v_i(ji,jj,jl) / a_i(ji,jj,jl) 
     461                     zhi = v_i(ji,jj,jl) / a_i(ji,jj,jl) 
    467462                     !---------------------------- 
    468463                     ! PE loss from deforming ice 
    469464                     !---------------------------- 
    470                      strength(ji,jj) = strength(ji,jj) - athorn(ji,jj,jl) * hi * hi 
     465                     strength(ji,jj) = strength(ji,jj) - athorn(ji,jj,jl) * zhi * zhi 
    471466 
    472467                     !-------------------------- 
    473468                     ! PE gain from rafting ice 
    474469                     !-------------------------- 
    475                      strength(ji,jj) = strength(ji,jj) + 2._wp * araft(ji,jj,jl) * hi * hi 
     470                     strength(ji,jj) = strength(ji,jj) + 2._wp * araft(ji,jj,jl) * zhi * zhi 
    476471 
    477472                     !---------------------------- 
     
    479474                     !---------------------------- 
    480475                     strength(ji,jj) = strength(ji,jj) + aridge(ji,jj,jl)/krdg(ji,jj,jl)     & 
    481                         * z1_3 * (hrmax(ji,jj,jl)**3 - hrmin(ji,jj,jl)**3) / ( hrmax(ji,jj,jl)-hrmin(ji,jj,jl) )    
    482 !!gm Optimization:  (a**3-b**3)/(a-b) = a*a+ab+b*b   ==> less costly operations even if a**3 is replaced by a*a*a...                     
     476                        * z1_3 * ( hrmax(ji,jj,jl)**2 + hrmin(ji,jj,jl)**2 +  hrmax(ji,jj,jl) * hrmin(ji,jj,jl) )   
     477                        !!(a**3-b**3)/(a-b) = a*a+ab+b*b                       
    483478                  ENDIF 
    484479                  ! 
     
    486481            END DO 
    487482         END DO 
    488  
    489          zzc = rn_pe_rdg * Cp     ! where Cp = (g/2)*(rhow-rhoi)*(rhoi/rhow) and rn_pe_rdg accounts for frictional dissipation 
    490          strength(:,:) = zzc * strength(:,:) / aksum(:,:) 
    491  
     483    
     484         strength(:,:) = rn_pe_rdg * Cp * strength(:,:) / aksum(:,:) 
     485                         ! where Cp = (g/2)*(rhow-rhoi)*(rhoi/rhow) and rn_pe_rdg accounts for frictional dissipation 
    492486         ksmooth = 1 
    493487 
     
    513507         DO jj = 1, jpj 
    514508            DO ji = 1, jpi 
    515                IF ( bv_i(ji,jj) > 0.0 ) THEN 
    516                   zdummy = MIN ( bv_i(ji,jj), 0.10 ) * MIN( bv_i(ji,jj), 0.10 ) 
    517                ELSE 
    518                   zdummy = 0.0 
    519                ENDIF 
    520509               strength(ji,jj) = strength(ji,jj) * exp(-5.88*SQRT(MAX(bv_i(ji,jj),0.0))) 
    521510            END DO 
     
    536525         CALL lbc_lnk( strength, 'T', 1. ) 
    537526 
    538          DO jj = 2, jpj - 1 
    539             DO ji = 2, jpi - 1 
     527         DO jj = 2, jpjm1 
     528            DO ji = 2, jpim1 
    540529               IF ( ( asum(ji,jj) - ato_i(ji,jj) ) > epsi10) THEN ! ice is present 
    541530                  zworka(ji,jj) = 4.0 * strength(ji,jj)              & 
     
    545534                     &          + strength(ji,jj+1) * tmask(ji,jj+1,1)     
    546535 
    547                   zw1 = 4.0 + tmask(ji-1,jj,1) + tmask(ji+1,jj,1) + tmask(ji,jj-1,1) + tmask(ji,jj+1,1) 
    548                   zworka(ji,jj) = zworka(ji,jj) / zw1 
     536                  zworka(ji,jj) = zworka(ji,jj) /  & 
     537                     &           ( 4.0 + tmask(ji-1,jj,1) + tmask(ji+1,jj,1) + tmask(ji,jj-1,1) + tmask(ji,jj+1,1) ) 
    549538               ELSE 
    550539                  zworka(ji,jj) = 0._wp 
     
    553542         END DO 
    554543 
    555          DO jj = 2, jpj - 1 
    556             DO ji = 2, jpi - 1 
     544         DO jj = 2, jpjm1 
     545            DO ji = 2, jpim1 
    557546               strength(ji,jj) = zworka(ji,jj) 
    558547            END DO 
     
    609598      !!---------------------------------------------------------------------! 
    610599      INTEGER ::   ji,jj, jl    ! dummy loop indices 
    611       REAL(wp) ::   Gstari, astari, hi, hrmean, zdummy   ! local scalar 
     600      REAL(wp) ::   Gstari, astari, zhi, hrmean, zdummy   ! local scalar 
    612601      REAL(wp), POINTER, DIMENSION(:,:)   ::   zworka    ! temporary array used here 
    613602      REAL(wp), POINTER, DIMENSION(:,:,:) ::   Gsum      ! Gsum(n) = sum of areas in categories 0 to n 
     
    697686                     athorn(ji,jj,jl) = 0.0 
    698687                  ENDIF 
    699                END DO ! ji 
    700             END DO ! jj 
    701          END DO ! jl  
     688               END DO 
     689            END DO 
     690         END DO 
    702691 
    703692      ELSE                             !--- Exponential, more stable formulation (Lipscomb et al, 2007) 
     
    791780 
    792781               IF (a_i(ji,jj,jl) > epsi10 .AND. athorn(ji,jj,jl) > 0.0 ) THEN 
    793                   hi = v_i(ji,jj,jl) / a_i(ji,jj,jl) 
    794                   hrmean          = MAX(SQRT(rn_hstar*hi), hi*krdgmin) 
    795                   hrmin(ji,jj,jl) = MIN(2.0*hi, 0.5*(hrmean + hi)) 
     782                  zhi = v_i(ji,jj,jl) / a_i(ji,jj,jl) 
     783                  hrmean          = MAX(SQRT(rn_hstar*zhi), zhi*krdgmin) 
     784                  hrmin(ji,jj,jl) = MIN(2.0*zhi, 0.5*(hrmean + zhi)) 
    796785                  hrmax(ji,jj,jl) = 2.0*hrmean - hrmin(ji,jj,jl) 
    797                   hraft(ji,jj,jl) = kraft*hi 
    798                   krdg(ji,jj,jl)  = hrmean / hi 
     786                  hraft(ji,jj,jl) = kraft*zhi 
     787                  krdg(ji,jj,jl)  = hrmean / zhi 
    799788               ELSE 
    800789                  hraft(ji,jj,jl) = 0.0 
     
    844833      INTEGER ::   ij                ! horizontal index, combines i and j loops 
    845834      INTEGER ::   icells            ! number of cells with aicen > puny 
    846       REAL(wp) ::   hL, hR, farea, zdummy, zdummy0, ztmelts    ! left and right limits of integration 
     835      REAL(wp) ::   hL, hR, farea, ztmelts    ! left and right limits of integration 
    847836      REAL(wp) ::   zsstK            ! SST in Kelvin 
    848837 
Note: See TracChangeset for help on using the changeset viewer.