Changeset 7421 for branches/2016/dev_merge_2016/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90

Timestamp:

2016-12-01T18:10:41+01:00 (8 years ago)

Author:

flavoni

Message:

#1811 merge dev_CNRS_MERATOR_2016 with dev_merge_2016 branch

File:

• branches/2016/dev_merge_2016/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90 (modified) (12 diffs)

branches/2016/dev_merge_2016/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90

@@ -21 +21 @@
    USE sbc_ice        ! Surface boundary condition: ice fields
    USE thd_ice        ! LIM thermodynamics
-   USE dom_ice        ! LIM domain
    USE ice            ! LIM variables
    USE limtab         ! LIM 2D <==> 1D
@@ -71 +70 @@
       !!               update ht_s_1d, ht_i_1d and tbif_1d(:,:)      
       !!------------------------------------------------------------------------
-      INTEGER ::   ji,jj,jk,jl      ! dummy loop indices
-      INTEGER ::   nbpac            ! local integers 
-      INTEGER ::   ii, ij, iter     !   -       -
-      REAL(wp)  ::   ztmelts, zdv, zfrazb, zweight, zde                         ! local scalars
+      INTEGER  ::   ji,jj,jk,jl      ! dummy loop indices
+      INTEGER  ::   nbpac            ! local integers 
+      INTEGER  ::   ii, ij, iter     !   -       -
+      REAL(wp) ::   ztmelts, zdv, zfrazb, zweight, zde                          ! local scalars
       REAL(wp) ::   zgamafr, zvfrx, zvgx, ztaux, ztwogp, zf                     !   -      -
       REAL(wp) ::   ztenagm, zvfry, zvgy, ztauy, zvrel2, zfp, zsqcd , zhicrit   !   -      -
@@ -154 +153 @@
 
       ! Default new ice thickness
-      WHERE( qlead(:,:) < 0._wp ) ; hicol = rn_hnewice
-      ELSEWHERE                   ; hicol = 0._wp
+      WHERE( qlead(:,:) < 0._wp ) ; hicol(:,:) = rn_hnewice
+      ELSEWHERE                   ; hicol(:,:) = 0._wp
       END WHERE
 
@@ -170 +169 @@
          zgamafr = 0.03
 
-         DO jj = 2, jpj
-            DO ji = 2, jpi
-               IF ( qlead(ji,jj) < 0._wp ) THEN
+         DO jj = 2, jpjm1
+            DO ji = 2, jpim1
+               IF ( qlead(ji,jj) < 0._wp .AND. tau_icebfr(ji,jj) == 0._wp ) THEN ! activated if cooling and no landfast
                   !-------------
                   ! Wind stress
@@ -195 +194 @@
                   !-------------------
                   ! C-grid ice velocity
-                  rswitch = MAX(  0._wp, SIGN( 1._wp , at_i(ji,jj) )  )
-                  zvgx    = rswitch * ( u_ice(ji-1,jj  ) * umask(ji-1,jj  ,1)  + u_ice(ji,jj) * umask(ji,jj,1) ) * 0.5_wp
-                  zvgy    = rswitch * ( v_ice(ji  ,jj-1) * vmask(ji  ,jj-1,1)  + v_ice(ji,jj) * vmask(ji,jj,1) ) * 0.5_wp
+                  zvgx    = ( u_ice(ji-1,jj  ) * umask(ji-1,jj  ,1)  + u_ice(ji,jj) * umask(ji,jj,1) ) * 0.5_wp
+                  zvgy    = ( v_ice(ji  ,jj-1) * vmask(ji  ,jj-1,1)  + v_ice(ji,jj) * vmask(ji,jj,1) ) * 0.5_wp
 
                   !-----------------------------------
@@ -203 +201 @@
                   !-----------------------------------
                   ! absolute relative velocity
-                  zvrel2 = MAX(  ( zvfrx - zvgx ) * ( zvfrx - zvgx )   &
-                     &         + ( zvfry - zvgy ) * ( zvfry - zvgy ) , 0.15 * 0.15 )
+                  rswitch      = MAX( 0._wp, SIGN( 1._wp , at_i(ji,jj) - epsi10 ) )
+                  zvrel2       = MAX(  ( zvfrx - zvgx ) * ( zvfrx - zvgx )   &
+                     &               + ( zvfry - zvgy ) * ( zvfry - zvgy ) , 0.15 * 0.15 ) * rswitch
                   zvrel(ji,jj) = SQRT( zvrel2 )
 
@@ -219 +218 @@
                      zfp = ( hicol(ji,jj) - zhicrit ) * ( 3.0 * hicol(ji,jj) + zhicrit ) - zhicrit * ztwogp * zvrel2
 
-                     hicol(ji,jj) = hicol(ji,jj) - zf/zfp
+                     hicol(ji,jj) = hicol(ji,jj) - zf / MAX( zfp, epsi20 )
                      iter = iter + 1
                   END DO
@@ -228 +227 @@
          END DO 
          ! 
-         CALL lbc_lnk( zvrel(:,:), 'T', 1. )
-         CALL lbc_lnk( hicol(:,:), 'T', 1. )
+         CALL lbc_lnk( zvrel, 'T', 1. )
+         CALL lbc_lnk( hicol, 'T', 1. )
 
       ENDIF ! End of computation of frazil ice collection thickness
@@ -240 +239 @@
       ! Select points for new ice formation
       !-------------------------------------
-      ! This occurs if open water energy budget is negative
+      ! This occurs if open water energy budget is negative (cooling) and there is no landfast ice
       nbpac = 0
       npac(:) = 0
@@ -246 +245 @@
       DO jj = 1, jpj
          DO ji = 1, jpi
-            IF ( qlead(ji,jj)  <  0._wp ) THEN
+            IF ( qlead(ji,jj)  <  0._wp .AND. tau_icebfr(ji,jj) == 0._wp ) THEN
                nbpac = nbpac + 1
                npac( nbpac ) = (jj - 1) * jpi + ji
@@ -255 +254 @@
       ! debug point to follow
       jiindex_1d = 0
-      IF( ln_icectl ) THEN
+      IF( ln_limctl ) THEN
          DO ji = mi0(iiceprt), mi1(iiceprt)
             DO jj = mj0(jiceprt), mj1(jiceprt)
@@ -265 +264 @@
       ENDIF
    
-      IF( ln_icectl ) WRITE(numout,*) 'lim_thd_lac : nbpac = ', nbpac
+      IF( ln_limctl ) WRITE(numout,*) 'lim_thd_lac : nbpac = ', nbpac
 
       !------------------------------