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 8813 for branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90 – NEMO

Ignore:
Timestamp:
2017-11-24T17:56:51+01:00 (6 years ago)
Author:
gm
Message:

#1911 (ENHANCE-09): PART I.3 - phasing with updated branch dev_r8183_ICEMODEL revision 8787

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r7881_ENHANCE09_RK3/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90

    r8637 r8813  
    4444   !!   ice_var_salprof1d : salinity profile in the ice 1D 
    4545   !!   ice_var_zapsmall  : remove very small area and volume 
    46    !!   ice_var_itd       : convert 1-cat to multiple cat 
     46   !!   ice_var_itd       : convert 1-cat to jpl-cat 
     47   !!   ice_var_itd2      : convert N-cat to jpl-cat 
    4748   !!   ice_var_bv        : brine volume 
    4849   !!---------------------------------------------------------------------- 
     
    6768   PUBLIC   ice_var_zapsmall 
    6869   PUBLIC   ice_var_itd 
     70   PUBLIC   ice_var_itd2 
    6971   PUBLIC   ice_var_bv            
    7072 
     
    9698      et_s(:,:)  = SUM( SUM( e_s(:,:,:,:), dim=4 ), dim=3 ) 
    9799      et_i(:,:)  = SUM( SUM( e_i(:,:,:,:), dim=4 ), dim=3 ) 
    98  
     100      ! 
    99101      at_ip(:,:) = SUM( a_ip(:,:,:), dim=3 ) ! melt ponds 
    100102      vt_ip(:,:) = SUM( v_ip(:,:,:), dim=3 ) 
    101  
    102       ato_i(:,:) = 1._wp - at_i(:,:)    ! open water fraction   
     103      ! 
     104      ato_i(:,:) = 1._wp - at_i(:,:)         ! open water fraction   
    103105 
    104106      IF( kn > 1 ) THEN 
     
    238240         END WHERE 
    239241      END DO 
    240  
     242      ! 
    241243      ! integrated values  
    242244      vt_i (:,:) = SUM( v_i, dim=3 ) 
     
    322324            END DO 
    323325         END DO 
    324  
     326         ! 
    325327         ! Computation of the profile 
    326328         DO jl = 1, jpl 
     
    362364   END SUBROUTINE ice_var_salprof 
    363365 
     366 
    364367   SUBROUTINE ice_var_salprof1d 
    365368      !!------------------------------------------------------------------- 
     
    457460         ELSEWHERE                                                                           ;   zswitch(:,:) = 1._wp 
    458461         END WHERE 
    459           
     462         ! 
    460463         DO jk = 1, nlay_i 
    461464            DO jj = 1 , jpj 
     
    468471            END DO 
    469472         END DO 
    470  
     473         ! 
    471474         DO jj = 1 , jpj 
    472475            DO ji = 1 , jpi 
     
    484487               t_s(ji,jj,1,jl) = t_s(ji,jj,1,jl) * zswitch(ji,jj) + rt0 * ( 1._wp - zswitch(ji,jj) ) 
    485488               e_s(ji,jj,1,jl) = e_s(ji,jj,1,jl) * zswitch(ji,jj) 
    486  
     489               ! 
    487490               !----------------------------------------------------------------- 
    488491               ! zap ice and snow volume, add water and salt to ocean 
     
    495498               oa_i (ji,jj,jl) = oa_i(ji,jj,jl) * zswitch(ji,jj) 
    496499               sv_i (ji,jj,jl) = sv_i(ji,jj,jl) * zswitch(ji,jj) 
    497  
     500               ! 
    498501               h_i (ji,jj,jl) = h_i (ji,jj,jl) * zswitch(ji,jj) 
    499502               h_s (ji,jj,jl) = h_s (ji,jj,jl) * zswitch(ji,jj) 
    500  
     503               ! 
    501504               a_ip (ji,jj,jl) = a_ip (ji,jj,jl) * zswitch(ji,jj) 
    502505               v_ip (ji,jj,jl) = v_ip (ji,jj,jl) * zswitch(ji,jj) 
    503  
    504             END DO 
    505          END DO 
    506          
     506               ! 
     507            END DO 
     508         END DO 
     509         ! 
    507510      END DO  
    508511 
     
    548551      !!------------------------------------------------------------------- 
    549552      INTEGER  :: ji, jk, jl             ! dummy loop indices 
    550       INTEGER  :: ijpij, i_fill, jl0   
     553      INTEGER  :: idim, i_fill, jl0   
    551554      REAL(wp) :: zarg, zV, zconv, zdh, zdv 
    552555      REAL(wp), DIMENSION(:),   INTENT(in)    ::   zhti, zhts, zati    ! input ice/snow variables 
     
    561564      ! then we check whether the distribution fullfills 
    562565      ! volume and area conservation, positivity and ice categories bounds 
    563       ijpij = SIZE( zhti , 1 ) 
    564       zh_i(1:ijpij,1:jpl) = 0._wp 
    565       zh_s(1:ijpij,1:jpl) = 0._wp 
    566       za_i (1:ijpij,1:jpl) = 0._wp 
    567  
    568       DO ji = 1, ijpij 
    569           
     566      idim = SIZE( zhti , 1 ) 
     567      zh_i(1:idim,1:jpl) = 0._wp 
     568      zh_s(1:idim,1:jpl) = 0._wp 
     569      za_i(1:idim,1:jpl) = 0._wp 
     570      ! 
     571      DO ji = 1, idim 
     572         ! 
    570573         IF( zhti(ji) > 0._wp ) THEN 
    571  
     574            ! 
    572575            ! find which category (jl0) the input ice thickness falls into 
    573576            jl0 = jpl 
     
    578581               ENDIF 
    579582            END DO 
    580  
     583            ! 
    581584            itest(:) = 0 
    582585            i_fill   = jpl + 1                                            !------------------------------------ 
     
    586589               ! 
    587590               zh_i(ji,1:jpl) = 0._wp 
    588                za_i (ji,1:jpl) = 0._wp 
    589                itest(:)        = 0       
    590                 
     591               za_i(ji,1:jpl) = 0._wp 
     592               itest(:)       = 0       
     593               ! 
    591594               IF ( i_fill == 1 ) THEN      !-- case very thin ice: fill only category 1 
    592595                  zh_i(ji,1) = zhti(ji) 
     
    597600                     zh_i(ji,jl) = hi_mean(jl) 
    598601                  END DO 
    599                    
     602                  ! 
    600603                  ! concentration 
    601604                  za_i(ji,jl0) = zati(ji) / SQRT(REAL(jpl)) 
     
    606609                     ENDIF 
    607610                  END DO 
    608                    
     611                  ! 
    609612                  ! last category 
    610613                  za_i(ji,i_fill) = zati(ji) - SUM( za_i(ji,1:i_fill-1) ) 
    611614                  zV = SUM( za_i(ji,1:i_fill-1) * zh_i(ji,1:i_fill-1) ) 
    612615                  zh_i(ji,i_fill) = ( zhti(ji) * zati(ji) - zV ) / MAX( za_i(ji,i_fill), epsi10 )  
    613                    
     616                  ! 
    614617                  ! clem: correction if concentration of upper cat is greater than lower cat 
    615618                  !       (it should be a gaussian around jl0 but sometimes it is not) 
     
    622625                           za_i (ji,1:jl-1) = za_i(ji,1:jl-1) + zdv / MAX( REAL(jl-1) * zhti(ji), epsi10 ) 
    623626                        END IF 
    624                      ENDDO 
     627                     END DO 
    625628                  ENDIF 
    626                 
     629                  ! 
    627630               ENDIF 
    628              
     631               ! 
    629632               ! Compatibility tests 
    630633               zconv = ABS( zati(ji) - SUM( za_i(ji,1:jpl) ) )  
    631                IF ( zconv < epsi06 ) itest(1) = 1                                        ! Test 1: area conservation 
    632              
     634               IF ( zconv < epsi06 )   itest(1) = 1                                        ! Test 1: area conservation 
     635               ! 
    633636               zconv = ABS( zhti(ji)*zati(ji) - SUM( za_i(ji,1:jpl)*zh_i(ji,1:jpl) ) ) 
    634                IF ( zconv < epsi06 ) itest(2) = 1                                        ! Test 2: volume conservation 
    635                 
    636                IF ( zh_i(ji,i_fill) >= hi_max(i_fill-1) ) itest(3) = 1                  ! Test 3: thickness of the last category is in-bounds ? 
    637                 
     637               IF ( zconv < epsi06 )   itest(2) = 1                                        ! Test 2: volume conservation 
     638               ! 
     639               IF ( zh_i(ji,i_fill) >= hi_max(i_fill-1) )   itest(3) = 1                  ! Test 3: thickness of the last category is in-bounds ? 
     640               ! 
    638641               itest(4) = 1 
    639642               DO jl = 1, i_fill 
     
    642645               !                                         !---------------------------- 
    643646            END DO                                       ! end iteration on categories 
    644                !                                         !---------------------------- 
     647            !                                            !---------------------------- 
    645648         ENDIF 
    646649      END DO 
     
    648651      ! Add Snow in each category where za_i is not 0 
    649652      DO jl = 1, jpl 
    650          DO ji = 1, ijpij 
     653         DO ji = 1, idim 
    651654            IF( za_i(ji,jl) > 0._wp ) THEN 
    652655               zh_s(ji,jl) = zh_i(ji,jl) * ( zhts(ji) / zhti(ji) ) 
     
    661664      END DO 
    662665      ! 
    663     END SUBROUTINE ice_var_itd 
    664  
    665  
    666     SUBROUTINE ice_var_bv 
     666   END SUBROUTINE ice_var_itd 
     667 
     668 
     669   SUBROUTINE ice_var_itd2( zhti, zhts, zati, zh_i, zh_s, za_i ) 
     670      !!------------------------------------------------------------------- 
     671      !!                ***  ROUTINE ice_var_itd2   *** 
     672      !! 
     673      !! ** Purpose :  converting N-cat ice to jpl ice categories 
     674      !! 
     675      !!                  ice thickness distribution follows a gaussian law 
     676      !!               around the concentration of the most likely ice thickness 
     677      !!                           (similar as iceistate.F90) 
     678      !! 
     679      !! ** Method:   Iterative procedure 
     680      !!                 
     681      !!               1) Fill ice cat that correspond to input thicknesses 
     682      !!                  Find the lowest(jlmin) and highest(jlmax) cat that are filled 
     683      !! 
     684      !!               2) Expand the filling to the cat jlmin-1 and jlmax+1 
     685      !!                   by removing 25% ice area from jlmin and jlmax (resp.)  
     686      !!               
     687      !!               3) Expand the filling to the empty cat between jlmin and jlmax  
     688      !!                   by a) removing 25% ice area from the lower cat (ascendant loop jlmin=>jlmax) 
     689      !!                      b) removing 25% ice area from the higher cat (descendant loop jlmax=>jlmin) 
     690      !! 
     691      !! ** Arguments : zhti: N-cat ice thickness 
     692      !!                zhts: N-cat snow depth 
     693      !!                zati: N-cat ice concentration 
     694      !! 
     695      !! ** Output    : jpl-cat  
     696      !! 
     697      !!  (Example of application: BDY forcings when inputs have N-cat /= jpl)   
     698      !!------------------------------------------------------------------- 
     699      INTEGER  ::   ji, jl, jl1, jl2             ! dummy loop indices 
     700      INTEGER  ::   idim, icat   
     701      INTEGER, PARAMETER ::   ztrans = 0.25_wp 
     702      REAL(wp), DIMENSION(:,:), INTENT(in)    ::   zhti, zhts, zati    ! input ice/snow variables 
     703      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zh_i, zh_s, za_i    ! output ice/snow variables 
     704      INTEGER , DIMENSION(:,:), ALLOCATABLE   ::   jlfil, jlfil2 
     705      INTEGER , DIMENSION(:)  , ALLOCATABLE   ::   jlmax, jlmin 
     706      !!------------------------------------------------------------------- 
     707      ! 
     708      idim = SIZE( zhti, 1 ) 
     709      icat = SIZE( zhti, 2 ) 
     710      ! 
     711      ALLOCATE( jlfil(idim,jpl), jlfil2(idim,jpl) )       ! allocate arrays 
     712      ALLOCATE( jlmin(idim), jlmax(idim) ) 
     713 
     714      ! --- initialize output fields to 0 --- ! 
     715      zh_i(1:idim,1:jpl) = 0._wp 
     716      zh_s(1:idim,1:jpl) = 0._wp 
     717      za_i(1:idim,1:jpl) = 0._wp 
     718      ! 
     719      ! --- fill the categories --- ! 
     720      !     find where cat-input = cat-output and fill cat-output fields   
     721      jlmax(:) = 0 
     722      jlmin(:) = 999 
     723      jlfil(:,:) = 0 
     724      DO jl1 = 1, jpl 
     725         DO jl2 = 1, icat 
     726            DO ji = 1, idim 
     727               IF( hi_max(jl1-1) <= zhti(ji,jl2) .AND. hi_max(jl1) > zhti(ji,jl2) ) THEN 
     728                  ! fill the right category 
     729                  zh_i(ji,jl1) = zhti(ji,jl2) 
     730                  zh_s(ji,jl1) = zhts(ji,jl2) 
     731                  za_i(ji,jl1) = zati(ji,jl2) 
     732                  ! record categories that are filled 
     733                  jlmax(ji) = MAX( jlmax(ji), jl1 ) 
     734                  jlmin(ji) = MIN( jlmin(ji), jl1 ) 
     735                  jlfil(ji,jl1) = jl1 
     736               ENDIF 
     737            END DO 
     738         END DO 
     739      END DO 
     740      ! 
     741      ! --- fill the gaps between categories --- !   
     742      !     transfer from categories filled at the previous step to the empty ones in between 
     743      DO ji = 1, idim 
     744         jl1 = jlmin(ji) 
     745         jl2 = jlmax(ji) 
     746         IF( jl1 > 1 ) THEN 
     747            ! fill the lower cat (jl1-1) 
     748            za_i(ji,jl1-1) = ztrans * za_i(ji,jl1) 
     749            zh_i(ji,jl1-1) = hi_mean(jl1-1) 
     750            ! remove from cat jl1 
     751            za_i(ji,jl1  ) = ( 1._wp - ztrans ) * za_i(ji,jl1) 
     752         ENDIF 
     753         IF( jl2 < jpl ) THEN 
     754            ! fill the upper cat (jl2+1) 
     755            za_i(ji,jl2+1) = ztrans * za_i(ji,jl2) 
     756            zh_i(ji,jl2+1) = hi_mean(jl2+1) 
     757            ! remove from cat jl2 
     758            za_i(ji,jl2  ) = ( 1._wp - ztrans ) * za_i(ji,jl2) 
     759         ENDIF 
     760      END DO 
     761      ! 
     762      jlfil2(:,:) = jlfil(:,:)  
     763      ! fill categories from low to high 
     764      DO jl = 2, jpl-1 
     765         DO ji = 1, idim 
     766            IF( jlfil(ji,jl-1) /= 0 .AND. jlfil(ji,jl) == 0 ) THEN 
     767               ! fill high 
     768               za_i(ji,jl) = ztrans * za_i(ji,jl-1) 
     769               zh_i(ji,jl) = hi_mean(jl) 
     770               jlfil(ji,jl) = jl 
     771               ! remove low 
     772               za_i(ji,jl-1) = ( 1._wp - ztrans ) * za_i(ji,jl-1) 
     773            ENDIF 
     774         END DO 
     775      END DO 
     776      ! 
     777      ! fill categories from high to low 
     778      DO jl = jpl-1, 2, -1 
     779         DO ji = 1, idim 
     780            IF( jlfil2(ji,jl+1) /= 0 .AND. jlfil2(ji,jl) == 0 ) THEN 
     781               ! fill low 
     782               za_i(ji,jl) = za_i(ji,jl) + ztrans * za_i(ji,jl+1) 
     783               zh_i(ji,jl) = hi_mean(jl)  
     784               jlfil2(ji,jl) = jl 
     785               ! remove high 
     786               za_i(ji,jl+1) = ( 1._wp - ztrans ) * za_i(ji,jl+1) 
     787            ENDIF 
     788         END DO 
     789      END DO 
     790      ! 
     791      DEALLOCATE( jlfil, jlfil2 )      ! deallocate arrays 
     792      DEALLOCATE( jlmin, jlmax ) 
     793      ! 
     794   END SUBROUTINE ice_var_itd2 
     795 
     796 
     797   SUBROUTINE ice_var_bv 
    667798      !!------------------------------------------------------------------- 
    668799      !!                ***  ROUTINE ice_var_bv *** 
Note: See TracChangeset for help on using the changeset viewer.