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 8500 – NEMO

Changeset 8500


Ignore:
Timestamp:
2017-09-06T12:01:34+02:00 (6 years ago)
Author:
clem
Message:

changes in style - part3 - move output into proper files and correct a bug in debug mode

Location:
branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3
Files:
9 edited

Legend:

Unmodified
Added
Removed

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/ice.F90

    r8498 r8500  
    475475   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:)   ::   t_si          !: Temperature at Snow-ice interface (K)  
    476476   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   tm_si         !: mean temperature at the snow-ice interface (K)  
    477    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   diag_dmi_dyn  !: Change in ice mass due to ice dynamics (kg/m2/s) 
    478    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   diag_dms_dyn  !: Change in snow mass due to ice dynamics (kg/m2/s) 
    479477   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   diag_fc_bo    !: Bottom conduction flux (W/m2) 
    480478   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   diag_fc_su    !: Surface conduction flux (W/m2) 
    481    REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:)     ::   diag_shear    !: Maximum shear of sea-ice velocity field 
    482479 
    483480   ! 
    … …  
    597594      ii = ii + 1 
    598595      ALLOCATE( t_si (jpi,jpj,jpl)    , tm_si(jpi,jpj)        ,    &  
    599                 diag_dmi_dyn(jpi,jpj) , diag_dms_dyn(jpi,jpj) ,    & 
    600596                diag_fc_bo(jpi,jpj)   , diag_fc_su(jpi,jpj)   ,    & 
    601                 diag_shear(jpi,jpj)   , STAT = ierr(ii) ) 
     597                STAT = ierr(ii) ) 
    602598 
    603599      ice_alloc = MAXVAL( ierr(:) ) 

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceadv.F90

    r8486 r8500  
    2929   USE lib_fortran    ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
    3030   USE timing         ! Timing 
     31   USE iom            ! 
    3132 
    3233   IMPLICIT NONE 
    … …  
    6869      ! 
    6970      REAL(wp) ::    zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b 
    70       REAL(wp) ::    zdv, zda 
     71      REAL(wp) ::    zdv 
    7172      REAL(wp), DIMENSION(jpi,jpj)           ::   zatold, zeiold, zesold, zsmvold  
    7273      REAL(wp), DIMENSION(jpi,jpj,jpl)       ::   zhimax, zviold, zvsold 
    … …  
    8081      ! MV MP 2016 
    8182      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:)   ::   z0ap , z0vp 
    82       REAL(wp) ::   za_old 
    8383      ! END MV MP 2016 
    8484      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:,:) ::   z0ei 
    … …  
    111111      !--- Thickness correction init. --- ! 
    112112      zatold(:,:) = at_i 
    113       DO jl = 1, jpl 
    114          DO jj = 1, jpj 
    115             DO ji = 1, jpi 
    116                rswitch          = MAX( 0._wp , SIGN( 1._wp, a_i(ji,jj,jl) - epsi20 ) ) 
    117                ht_i  (ji,jj,jl) = v_i (ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi20 ) * rswitch 
    118                ht_s  (ji,jj,jl) = v_s (ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi20 ) * rswitch 
    119             END DO 
    120          END DO 
    121       END DO 
     113      WHERE( a_i(:,:,:) >= epsi20 ) 
     114         ht_i(:,:,:) = v_i(:,:,:) / a_i(:,:,:) 
     115         ht_s(:,:,:) = v_s(:,:,:) / a_i(:,:,:) 
     116      ELSEWHERE 
     117         ht_i(:,:,:) = 0._wp 
     118         ht_s(:,:,:) = 0._wp          
     119      END WHERE 
     120          
    122121      ! --- Record max of the surrounding ice thicknesses for correction in case advection creates ice too thick --- ! 
    123122      zhimax(:,:,:) = ht_i(:,:,:) + ht_s(:,:,:) 
    … …  
    126125            DO ji = 2, jpim1 
    127126!!gm use of MAXVAL here is very probably less efficient than expending the 9 values 
    128                zhimax(ji,jj,jl) = MAXVAL( ht_i(ji-1:ji+1,jj-1:jj+1,jl) + ht_s(ji-1:ji+1,jj-1:jj+1,jl) ) 
     127               zhimax(ji,jj,jl) = MAX( epsi20, MAXVAL( ht_i(ji-1:ji+1,jj-1:jj+1,jl) + ht_s(ji-1:ji+1,jj-1:jj+1,jl) ) ) 
    129128            END DO 
    130129         END DO 
    … …  
    399398         END DO 
    400399      END DO 
    401        
    402       IF( nn_limdyn == 2) THEN 
     400      IF( iom_use('icetrp') )   CALL iom_put( "icetrp" , diag_trp_vi * rday  )         ! ice volume transport 
     401      IF( iom_use('snwtrp') )   CALL iom_put( "snwtrp" , diag_trp_vs * rday  )         ! snw volume transport 
     402      IF( iom_use('saltrp') )   CALL iom_put( "saltrp" , diag_trp_smv * rday * rhoic ) ! salt content transport 
     403      IF( iom_use('deitrp') )   CALL iom_put( "deitrp" , diag_trp_ei         )         ! advected ice enthalpy (W/m2) 
     404      IF( iom_use('destrp') )   CALL iom_put( "destrp" , diag_trp_es         )         ! advected snw enthalpy (W/m2) 
     405       
     406      IF( nn_limdyn == 2 ) THEN 
    403407         ! 
    404408         CALL ice_var_zapsmall      !--- zap small areas ---! 
    … …  
    407411            DO jj = 1, jpj 
    408412               DO ji = 1, jpi 
    409                   !  
    410413                  IF ( v_i(ji,jj,jl) > 0._wp ) THEN 
    411414                     ! 
    412                      rswitch          = MAX( 0._wp , SIGN( 1._wp, a_i(ji,jj,jl) - epsi20 ) ) 
    413                      ht_i  (ji,jj,jl) = v_i (ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi20 ) * rswitch 
    414                      ht_s  (ji,jj,jl) = v_s (ji,jj,jl) / MAX( a_i(ji,jj,jl) , epsi20 ) * rswitch 
    415                      ! 
    416                      zdv  = v_i(ji,jj,jl) + v_s(ji,jj,jl) - zviold(ji,jj,jl) - zvsold(ji,jj,jl)   
     415                     ht_i  (ji,jj,jl) = v_i (ji,jj,jl) / a_i(ji,jj,jl) 
     416                     ht_s  (ji,jj,jl) = v_s (ji,jj,jl) / a_i(ji,jj,jl) 
     417                     zdv = v_i(ji,jj,jl) + v_s(ji,jj,jl) - zviold(ji,jj,jl) - zvsold(ji,jj,jl)   
    417418                     ! 
    418419                     IF ( ( zdv >  0.0 .AND. (ht_i(ji,jj,jl)+ht_s(ji,jj,jl)) > zhimax(ji,jj,jl) .AND. zatold(ji,jj) < 0.80 ) .OR. & 
    419420                        & ( zdv <= 0.0 .AND. (ht_i(ji,jj,jl)+ht_s(ji,jj,jl)) > zhimax(ji,jj,jl) ) ) THEN 
    420                         ! 
    421                         rswitch        = MAX( 0._wp, SIGN( 1._wp, zhimax(ji,jj,jl) - epsi20 ) ) 
    422                         a_i(ji,jj,jl)  = rswitch * ( v_i(ji,jj,jl) + v_s(ji,jj,jl) ) / MAX( zhimax(ji,jj,jl), epsi20 ) 
    423                         ! 
    424                         ! small correction due to *rswitch for a_i 
    425                         v_i  (ji,jj,jl)        = rswitch * v_i  (ji,jj,jl) 
    426                         v_s  (ji,jj,jl)        = rswitch * v_s  (ji,jj,jl) 
    427                         smv_i(ji,jj,jl)        = rswitch * smv_i(ji,jj,jl) 
    428                         e_s(ji,jj,1,jl)        = rswitch * e_s(ji,jj,1,jl) 
    429                         e_i(ji,jj,1:nlay_i,jl) = rswitch * e_i(ji,jj,1:nlay_i,jl) 
    430  
    431                         ! MV MP 2016 
    432                         IF ( nn_pnd_scheme > 0 ) THEN 
    433                            a_ip (ji,jj,jl)        = rswitch * a_ip (ji,jj,jl) 
    434                            v_ip (ji,jj,jl)        = rswitch * v_ip (ji,jj,jl) 
    435                         ENDIF 
    436                         ! END MV MP 2016 
    437                         ! 
     421                        a_i (ji,jj,jl) = ( v_i(ji,jj,jl) + v_s(ji,jj,jl) ) / zhimax(ji,jj,jl) 
     422                        ht_i(ji,jj,jl) =   v_i(ji,jj,jl) / a_i(ji,jj,jl) 
    438423                     ENDIF 
    439424                     ! 
    440425                  ENDIF 
    441                   ! 
    442426               END DO 
    443427            END DO 
    444428         END DO 
    445           
    446          DO jj = 1, jpj             !--- bound ht_i to hi_max (99 m). 
    447             DO ji = 1, jpi 
    448                ! MV MP 2016 
    449                za_old = a_i(ji,jj,jpl) 
    450                ! END MV MP 2016 
    451                rswitch         = MAX( 0._wp , SIGN( 1._wp, ht_i(ji,jj,jpl) - epsi20 ) ) 
    452                ht_i(ji,jj,jpl) = MIN( ht_i(ji,jj,jpl) , hi_max(jpl) ) 
    453                a_i (ji,jj,jpl) = v_i(ji,jj,jpl) / MAX( ht_i(ji,jj,jpl) , epsi20 ) * rswitch 
    454                ! MV MP 2016 
    455                IF ( nn_pnd_scheme > 0 ) THEN 
    456                   ! correct pond fraction to avoid a_ip > a_i 
    457                   a_ip(ji,jj,jpl) = a_ip(ji,jj,jpl) * a_i(ji,jj,jpl) / MAX( za_old , epsi20 ) * rswitch 
    458                ENDIF 
    459                ! END MP 2016 
    460             END DO 
    461          END DO 
     429                   
     430         WHERE( ht_i(:,:,jpl) > hi_max(jpl) )             !--- bound ht_i to hi_max (99 m) 
     431            ht_i(:,:,jpl) = hi_max(jpl) 
     432            a_i (:,:,jpl) = v_i(:,:,jpl) / hi_max(jpl) 
     433         END WHERE 
     434 
     435         IF ( nn_pnd_scheme > 0 ) THEN                    !--- correct pond fraction to avoid a_ip > a_i 
     436            WHERE( a_ip(:,:,:) > a_i(:,:,:) )   a_ip(:,:,:) = a_i(:,:,:) 
     437         ENDIF 
    462438         ! 
    463439      ENDIF 
    … …  
    467443      !------------------------------------------------------------ 
    468444      ! 
    469 !!gm remplace the test by, l_piling a logical compute one for all in icestp.F90 (and its declaration in ice.F90 
    470 !!gm      IF ( nn_limdyn == 1 .OR. ( ( nn_monocat == 2 ) .AND. ( jpl == 1 ) ) ) THEN ! simple conservative piling, comparable with LIM2 
    471       IF( l_piling ) THEN 
     445      IF( l_piling ) THEN ! simple conservative piling, comparable with 1-cat models 
    472446         at_i(:,:) = SUM( a_i(:,:,:), dim=3 ) 
    473447         DO jl = 1, jpl 
    474             DO jj = 1, jpj 
    475                DO ji = 1, jpi 
    476                   rswitch = MAX( 0._wp, SIGN( 1._wp, at_i(ji,jj) - epsi20 ) ) 
    477                   zda     = rswitch * MIN( rn_amax_2d(ji,jj) - at_i(ji,jj), 0._wp ) * a_i(ji,jj,jl) / MAX( at_i(ji,jj), epsi20 ) 
    478                   a_i(ji,jj,jl) = a_i(ji,jj,jl) + zda 
    479                END DO 
    480             END DO 
    481          END DO 
    482 !!gm better and faster coding? 
    483 !         DO jl = 1, jpl 
    484 !            WHERE( at_i(:,:) > epsi20 ) 
    485 !               a_i(:,:,jl) = a_i(:,:,jl) * (  1._wp + MIN( rn_amax_2d(:,:) - at_i(:,:) , 0._wp ) / at_i(:,:)  ) 
    486 !            END WHERE 
    487 !         END DO 
    488 !!gm end 
     448            WHERE( at_i(:,:) > epsi20 ) 
     449               a_i(:,:,jl) = a_i(:,:,jl) * (  1._wp + MIN( rn_amax_2d(:,:) - at_i(:,:) , 0._wp ) / at_i(:,:)  ) 
     450            END WHERE 
     451         END DO 
    489452      ENDIF 
    490453       

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/iceitd.F90

    r8498 r8500  
    422422            ! 
    423423            jl1 = kdonor(ji,jl) 
    424             IF ( jl1 == jl  ) THEN   ;   jl2 = jl1+1 
    425             ELSE                     ;   jl2 = jl  
    426             ENDIF 
    427             ! 
    428             IF( v_i_2d(ji,jl1) >= epsi10 ) THEN   ;   zworkv(ji) = pdvice(ji,jl) / v_i_2d(ji,jl1) 
    429             ELSE                                  ;   zworkv(ji) = 0._wp 
    430             ENDIF 
    431             IF( a_i_2d(ji,jl1) >= epsi10 ) THEN   ;   zworka(ji) = pdaice(ji,jl) / a_i_2d(ji,jl1) 
    432             ELSE                                  ;   zworka(ji) = 0._wp 
    433             ENDIF 
    434             ! 
    435             a_i_2d(ji,jl1) = a_i_2d(ji,jl1) - pdaice(ji,jl)       ! Ice areas 
    436             a_i_2d(ji,jl2) = a_i_2d(ji,jl2) + pdaice(ji,jl) 
    437             ! 
    438             v_i_2d(ji,jl1) = v_i_2d(ji,jl1) - pdvice(ji,jl)       ! Ice volumes 
    439             v_i_2d(ji,jl2) = v_i_2d(ji,jl2) + pdvice(ji,jl) 
    440             ! 
    441             ztrans         = v_s_2d(ji,jl1) * zworkv(ji)          ! Snow volumes 
    442             v_s_2d(ji,jl1) = v_s_2d(ji,jl1) - ztrans 
    443             v_s_2d(ji,jl2) = v_s_2d(ji,jl2) + ztrans  
    444             !           
    445             !                                                     ! Ice age  
    446             ztrans             = oa_i_2d(ji,jl1) * pdaice(ji,jl)  !!clem: should be * zworka(ji) but it does not work ???? 
    447             oa_i_2d(ji,jl1)    = oa_i_2d(ji,jl1) - ztrans 
    448             oa_i_2d(ji,jl2)    = oa_i_2d(ji,jl2) + ztrans 
    449             ! 
    450             ztrans             = smv_i_2d(ji,jl1) * zworkv(ji)    ! Ice salinity 
    451  
    452             smv_i_2d(ji,jl1)   = smv_i_2d(ji,jl1) - ztrans 
    453             smv_i_2d(ji,jl2)   = smv_i_2d(ji,jl2) + ztrans 
    454             ! 
    455             !                                                     ! Surface temperature 
    456             ztrans          = t_su_2d(ji,jl1) * pdaice(ji,jl)     !!clem: should be * zworka(ji) but it does not work ???? 
    457             zaTsfn(ji,jl1)  = zaTsfn(ji,jl1) - ztrans 
    458             zaTsfn(ji,jl2)  = zaTsfn(ji,jl2) + ztrans 
    459             !   
    460             ! MV MP 2016  
    461             IF ( nn_pnd_scheme > 0 ) THEN 
    462                !                                                  ! Pond fraction 
    463                ztrans          = a_ip_2d(ji,jl1) * pdaice(ji,jl) !!clem: should be * zworka(ji) but it does not work 
    464                a_ip_2d(ji,jl1) = a_ip_2d(ji,jl1) - ztrans 
    465                a_ip_2d(ji,jl2) = a_ip_2d(ji,jl2) + ztrans 
    466                !                                                  ! Pond volume (also proportional to da/a) 
    467                ztrans          = v_ip_2d(ji,jl1) * pdaice(ji,jl) !!clem: should be * zworka(ji) but it does not work 
    468                v_ip_2d(ji,jl1) = v_ip_2d(ji,jl1) - ztrans 
    469                v_ip_2d(ji,jl2) = v_ip_2d(ji,jl2) + ztrans 
    470             ENDIF 
    471             ! END MV MP 2016  
     424            ! 
     425            IF( jl1 > 0 ) THEN 
     426               ! 
     427               IF ( jl1 == jl  ) THEN   ;   jl2 = jl1+1 
     428               ELSE                     ;   jl2 = jl  
     429               ENDIF 
     430               ! 
     431               IF( v_i_2d(ji,jl1) >= epsi10 ) THEN   ;   zworkv(ji) = pdvice(ji,jl) / v_i_2d(ji,jl1) 
     432               ELSE                                  ;   zworkv(ji) = 0._wp 
     433               ENDIF 
     434               IF( a_i_2d(ji,jl1) >= epsi10 ) THEN   ;   zworka(ji) = pdaice(ji,jl) / a_i_2d(ji,jl1) 
     435               ELSE                                  ;   zworka(ji) = 0._wp 
     436               ENDIF 
     437               ! 
     438               a_i_2d(ji,jl1) = a_i_2d(ji,jl1) - pdaice(ji,jl)       ! Ice areas 
     439               a_i_2d(ji,jl2) = a_i_2d(ji,jl2) + pdaice(ji,jl) 
     440               ! 
     441               v_i_2d(ji,jl1) = v_i_2d(ji,jl1) - pdvice(ji,jl)       ! Ice volumes 
     442               v_i_2d(ji,jl2) = v_i_2d(ji,jl2) + pdvice(ji,jl) 
     443               ! 
     444               ztrans         = v_s_2d(ji,jl1) * zworkv(ji)          ! Snow volumes 
     445               v_s_2d(ji,jl1) = v_s_2d(ji,jl1) - ztrans 
     446               v_s_2d(ji,jl2) = v_s_2d(ji,jl2) + ztrans  
     447               !           
     448               !                                                     ! Ice age  
     449               ztrans             = oa_i_2d(ji,jl1) * pdaice(ji,jl)  !!clem: should be * zworka(ji) but it does not work ???? 
     450               oa_i_2d(ji,jl1)    = oa_i_2d(ji,jl1) - ztrans 
     451               oa_i_2d(ji,jl2)    = oa_i_2d(ji,jl2) + ztrans 
     452               ! 
     453               ztrans             = smv_i_2d(ji,jl1) * zworkv(ji)    ! Ice salinity 
     454               ! 
     455               smv_i_2d(ji,jl1)   = smv_i_2d(ji,jl1) - ztrans 
     456               smv_i_2d(ji,jl2)   = smv_i_2d(ji,jl2) + ztrans 
     457               ! 
     458               !                                                     ! Surface temperature 
     459               ztrans          = t_su_2d(ji,jl1) * pdaice(ji,jl)     !!clem: should be * zworka(ji) but it does not work ???? 
     460               zaTsfn(ji,jl1)  = zaTsfn(ji,jl1) - ztrans 
     461               zaTsfn(ji,jl2)  = zaTsfn(ji,jl2) + ztrans 
     462               !   
     463               ! MV MP 2016  
     464               IF ( nn_pnd_scheme > 0 ) THEN 
     465                  !                                                  ! Pond fraction 
     466                  ztrans          = a_ip_2d(ji,jl1) * pdaice(ji,jl) !!clem: should be * zworka(ji) but it does not work 
     467                  a_ip_2d(ji,jl1) = a_ip_2d(ji,jl1) - ztrans 
     468                  a_ip_2d(ji,jl2) = a_ip_2d(ji,jl2) + ztrans 
     469                  !                                                  ! Pond volume (also proportional to da/a) 
     470                  ztrans          = v_ip_2d(ji,jl1) * pdaice(ji,jl) !!clem: should be * zworka(ji) but it does not work 
     471                  v_ip_2d(ji,jl1) = v_ip_2d(ji,jl1) - ztrans 
     472                  v_ip_2d(ji,jl2) = v_ip_2d(ji,jl2) + ztrans 
     473               ENDIF 
     474               ! END MV MP 2016 
     475               ! 
     476            ENDIF   ! jl1 >0 
    472477         END DO 
    473478         ! 
    … …  
    479484               ! 
    480485               jl1 = kdonor(ji,jl) 
    481                IF(jl1 == jl) THEN  ;  jl2 = jl+1 
    482                ELSE                ;  jl2 = jl 
    483                ENDIF 
    484                ! 
    485                ztrans            = e_s(ii,ij,jk,jl1) * zworkv(ji) 
    486                e_s(ii,ij,jk,jl1) = e_s(ii,ij,jk,jl1) - ztrans 
    487                e_s(ii,ij,jk,jl2) = e_s(ii,ij,jk,jl2) + ztrans  
     486               ! 
     487               IF( jl1 > 0 ) THEN 
     488                  IF(jl1 == jl) THEN  ;  jl2 = jl+1 
     489                  ELSE                ;  jl2 = jl 
     490                  ENDIF 
     491                  ! 
     492                  ztrans            = e_s(ii,ij,jk,jl1) * zworkv(ji) 
     493                  e_s(ii,ij,jk,jl1) = e_s(ii,ij,jk,jl1) - ztrans 
     494                  e_s(ii,ij,jk,jl2) = e_s(ii,ij,jk,jl2) + ztrans 
     495               ENDIF 
    488496            END DO 
    489497         END DO 
    … …  
    495503               ! 
    496504               jl1 = kdonor(ji,jl) 
    497                IF(jl1 == jl) THEN  ;  jl2 = jl+1 
    498                ELSE                ;  jl2 = jl 
    499                ENDIF 
    500                ! 
    501                ztrans            = e_i(ii,ij,jk,jl1) * zworkv(ji) 
    502                e_i(ii,ij,jk,jl1) = e_i(ii,ij,jk,jl1) - ztrans 
    503                e_i(ii,ij,jk,jl2) = e_i(ii,ij,jk,jl2) + ztrans  
     505               ! 
     506               IF( jl1 > 0 ) THEN 
     507                  IF(jl1 == jl) THEN  ;  jl2 = jl+1 
     508                  ELSE                ;  jl2 = jl 
     509                  ENDIF 
     510                  ! 
     511                  ztrans            = e_i(ii,ij,jk,jl1) * zworkv(ji) 
     512                  e_i(ii,ij,jk,jl1) = e_i(ii,ij,jk,jl1) - ztrans 
     513                  e_i(ii,ij,jk,jl2) = e_i(ii,ij,jk,jl2) + ztrans 
     514               ENDIF 
    504515            END DO 
    505516         END DO 
    … …  
    510521      ! 3) Update ice thickness and temperature 
    511522      !------------------------------------------------------------------------------- 
    512       DO jl = 1, jpl 
    513          DO ji = 1, nidx 
    514             IF ( a_i_2d(ji,jl) > epsi10 ) THEN  
    515                ht_i_2d(ji,jl)  =  v_i_2d(ji,jl) / a_i_2d(ji,jl)  
    516                t_su_2d(ji,jl)  =  zaTsfn(ji,jl) / a_i_2d(ji,jl)  
    517             ELSE 
    518                ht_i_2d(ji,jl)  = 0._wp 
    519                t_su_2d(ji,jl)  = rt0 
    520             ENDIF 
    521          END DO 
    522       END DO 
     523      WHERE( a_i_2d(1:nidx,:) >= epsi20 ) 
     524         ht_i_2d(1:nidx,:)  =  v_i_2d(1:nidx,:) / a_i_2d(1:nidx,:)  
     525         t_su_2d(1:nidx,:)  =  zaTsfn(1:nidx,:) / a_i_2d(1:nidx,:)  
     526      ELSEWHERE 
     527         ht_i_2d(1:nidx,:)  = 0._wp 
     528         t_su_2d(1:nidx,:)  = rt0 
     529      END WHERE 
    523530      ! 
    524531      CALL tab_2d_3d( nidx, idxice(1:nidx), ht_i_2d (1:nidx,1:jpl), ht_i  ) 

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icerdgrft.F90

    r8498 r8500  
    612612               &                                      + rhosn * vsrft(ij) * ( 1._wp - rn_fsnowrft ) ) * r1_rdtice  ! fresh water source for ocean 
    613613 
    614 !!!            wfx_snw_dyn(ji,jj) = wfx_snw_dyn(ji,jj) + zwfx_snw 
    615  
    616614            hfx_dyn(ji,jj) = hfx_dyn(ji,jj) + ( - esrdg(ij) * ( 1._wp - rn_fsnowrdg )         &  
    617615               &                                - esrft(ij) * ( 1._wp - rn_fsnowrft ) ) * r1_rdtice        ! heat sink for ocean (<0, W.m-2) 
    … …  
    734732         ! 
    735733      END DO ! jl1 (deforming categories) 
    736  
    737       ! SIMIP diagnostics 
    738       diag_dmi_dyn(:,:) = - wfx_dyn(:,:)     + rhoic * diag_trp_vi(:,:) 
    739       diag_dms_dyn(:,:) = - wfx_snw_dyn(:,:) + rhosn * diag_trp_vs(:,:)       
    740734      ! 
    741735   END SUBROUTINE ice_rdgrft_ridgeshift 

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icerhg_evp.F90

    r8498 r8500  
    153153      REAL(wp), PARAMETER          ::   zmmin  = 1._wp                  ! ice mass (kg/m2) below which ice velocity equals ocean velocity 
    154154      !! --- diags 
    155       REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zswi, zsig1, zsig2, zsig3 
     155      REAL(wp), DIMENSION(jpi,jpj) ::   zswi 
     156      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zsig1, zsig2, zsig3 
    156157      !! --- SIMIP diags 
    157158      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   zdiag_sig1      ! Average normal stress in sea ice    
    … …  
    692693      ! 5) diagnostics 
    693694      !------------------------------------------------------------------------------! 
    694       ! --- ellipse --- ! 
     695      DO jj = 1, jpj 
     696         DO ji = 1, jpi 
     697            zswi(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice, 0 if no ice 
     698         END DO 
     699      END DO 
     700 
     701      ! --- divergence, shear and strength --- ! 
     702      IF( iom_use('idive')  )   CALL iom_put( "idive"  , divu_i(:,:)   * zswi(:,:) )   ! divergence 
     703      IF( iom_use('ishear') )   CALL iom_put( "ishear" , shear_i(:,:)  * zswi(:,:) )   ! shear 
     704      IF( iom_use('icestr') )   CALL iom_put( "icestr" , strength(:,:) * zswi(:,:) )   ! Ice strength 
     705 
     706      ! --- charge ellipse --- ! 
    695707      IF( iom_use('isig1') .OR. iom_use('isig2') .OR. iom_use('isig3') ) THEN 
    696708         ! 
    697          ALLOCATE( zswi(jpi,jpj) , zsig1(jpi,jpj) , zsig2(jpi,jpj) , zsig3(jpi,jpj) ) 
    698          ! 
    699          DO jj = 1, jpj 
    700             DO ji = 1, jpi 
    701                zswi(ji,jj) = MAX( 0._wp , SIGN( 1._wp , at_i(ji,jj) - epsi06 ) ) ! 1 if ice, 0 if no ice 
    702             END DO 
    703          END DO 
    704           
     709         ALLOCATE( zsig1(jpi,jpj) , zsig2(jpi,jpj) , zsig3(jpi,jpj) ) 
     710         !          
    705711         DO jj = 2, jpjm1 
    706712            DO ji = 2, jpim1 
    … …  
    728734         IF( iom_use('isig3') )   CALL iom_put( "isig3" , zsig3 ) 
    729735         ! 
    730          DEALLOCATE( zswi , zsig1 , zsig2 , zsig3 ) 
     736         DEALLOCATE( zsig1 , zsig2 , zsig3 ) 
    731737      ENDIF 
    732738       

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icestp.F90

    r8498 r8500  
    571571             
    572572            ! SIMIP diagnostics 
    573             diag_dms_dyn(ji,jj)  = 0._wp ; diag_dmi_dyn(ji,jj)  = 0._wp 
    574573            diag_fc_bo(ji,jj)    = 0._wp ; diag_fc_su(ji,jj)    = 0._wp 
    575574             

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icethd_dif.F90

    r8486 r8500  
    738738      ! MV SIMIP 2016 
    739739      !--- Conduction fluxes (positive downwards) 
    740       diag_fc_bo_1d(:)   = diag_fc_bo_1d(:) + fc_bo_i(:) * a_i_1d(:) / at_i_1d(:) 
    741       diag_fc_su_1d(:)   = diag_fc_su_1d(:) + fc_su(:)   * a_i_1d(:) / at_i_1d(:) 
     740      DO ji = 1, nidx 
     741         diag_fc_bo_1d(ji)   = diag_fc_bo_1d(ji) + fc_bo_i(ji) * a_i_1d(ji) / at_i_1d(ji) 
     742         diag_fc_su_1d(ji)   = diag_fc_su_1d(ji) + fc_su(ji)   * a_i_1d(ji) / at_i_1d(ji) 
     743      END DO 
    742744      ! END MV SIMIP 2016 
    743745 

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icevar.F90

    r8498 r8500  
    116116         ! 
    117117         ALLOCATE( z1_at_i(jpi,jpj) , z1_vt_i(jpi,jpj) ) 
    118          WHERE( at_i(:,:) > epsi10 )   ;   z1_at_i(:,:) = 1._wp / at_i(:,:) 
     118         WHERE( at_i(:,:) > epsi20 )   ;   z1_at_i(:,:) = 1._wp / at_i(:,:) 
    119119         ELSEWHERE                     ;   z1_at_i(:,:) = 0._wp 
    120120         END WHERE 
    121          WHERE( vt_i(:,:) > epsi10 )   ;   z1_vt_i(:,:) = 1._wp / vt_i(:,:) 
     121         WHERE( vt_i(:,:) > epsi20 )   ;   z1_vt_i(:,:) = 1._wp / vt_i(:,:) 
    122122         ELSEWHERE                     ;   z1_vt_i(:,:) = 0._wp 
    123123         END WHERE 
    … …  
    163163      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   z1_a_i, z1_v_i 
    164164      !!------------------------------------------------------------------ 
    165  
    166 !!gm Question 1:  here use epsi20 , in ice_var_agg it is epsi10 which is used...  why ?? 
    167165 
    168166!!gm Question 2:  It is possible to define existence of sea-ice in a common way between  
    … …  
    549547      ! to be sure that at_i is the sum of a_i(jl) 
    550548      at_i (:,:) = a_i(:,:,1) 
    551      vt_i (:,:) = v_i(:,:,1) 
     549      vt_i (:,:) = v_i(:,:,1) 
    552550      DO jl = 2, jpl 
    553551         at_i(:,:) = at_i(:,:) + a_i(:,:,jl) 
    554         vt_i(:,:) = vt_i(:,:) + v_i(:,:,jl) 
     552         vt_i(:,:) = vt_i(:,:) + v_i(:,:,jl) 
    555553      END DO 
    556554 

  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/icewri.F90

    r8498 r8500  
    145145      CALL iom_put( "micesalt"    , smt_i   * zswi      )        ! mean ice salinity 
    146146      CALL iom_put( "snowvol"     , vt_s    * zswi      )        ! snow volume 
    147  
    148       CALL iom_put( "idive"       , divu_i(:,:)  * zswi(:,:) + zmiss(:,:) )   ! divergence 
    149       CALL iom_put( "ishear"      , shear_i(:,:) * zswi(:,:) + zmiss(:,:) )   ! shear 
    150147       
    151       CALL iom_put( "icetrp"      , diag_trp_vi * rday  )        ! ice volume transport 
    152       CALL iom_put( "snwtrp"      , diag_trp_vs * rday  )        ! snw volume transport 
    153       CALL iom_put( "saltrp"      , diag_trp_smv * rday * rhoic ) ! salt content transport 
    154       CALL iom_put( "deitrp"      , diag_trp_ei         )        ! advected ice enthalpy (W/m2) 
    155       CALL iom_put( "destrp"      , diag_trp_es         )        ! advected snw enthalpy (W/m2) 
    156  
    157148      CALL iom_put( "sfxbog"      , sfx_bog * rday      )        ! salt flux from bottom growth 
    158149      CALL iom_put( "sfxbom"      , sfx_bom * rday      )        ! salt flux from bottom melting 
    … …  
    275266      IF ( iom_use( "dmithd"   ) ) CALL iom_put( "dmithd"      , - wfx_bog - wfx_bom - wfx_sum   &                       ! Sea-ice mass change from thermodynamics 
    276267              &                                                   - wfx_sni - wfx_opw - wfx_res ) 
    277       IF ( iom_use( "dmidyn"   ) ) CALL iom_put( "dmidyn"      ,   diag_dmi_dyn             )                            ! Sea-ice mass change from dynamics 
     268      IF ( iom_use( "dmidyn"   ) ) CALL iom_put( "dmidyn"      , - wfx_dyn(:,:) + rhoic * diag_trp_vi(:,:) )             ! Sea-ice mass change from dynamics(kg/m2/s) 
    278269      IF ( iom_use( "dmiopw"   ) ) CALL iom_put( "dmiopw"      , - wfx_opw                  )                            ! Sea-ice mass change through growth in open water 
    279270      IF ( iom_use( "dmibog"   ) ) CALL iom_put( "dmibog"      , - wfx_bog                  )                            ! Sea-ice mass change through basal growth 
    … …  
    290281 
    291282      IF ( iom_use( "dmsmel"   ) ) CALL iom_put( "dmsmel"      , - wfx_snw_sum              )                            ! Snow mass change through melt 
    292       IF ( iom_use( "dmsdyn"   ) ) CALL iom_put( "dmsdyn"      ,   diag_dms_dyn             )                            ! Snow mass change through dynamics 
     283      IF ( iom_use( "dmsdyn"   ) ) CALL iom_put( "dmsdyn"      , - wfx_snw_dyn(:,:) + rhosn * diag_trp_vs(:,:) )         ! Snow mass change through dynamics(kg/m2/s) 
    293284 
    294285      IF ( iom_use( "hfxsenso" ) ) CALL iom_put( "hfxsenso"    ,   -fhtur(:,:)              * zswi(:,:) + zmiss(:,:) )   ! Sensible oceanic heat flux 
    … …  
    305296      IF ( iom_use( "utau_ice" ) ) CALL iom_put( "utau_ice"     ,  utau_ice(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Wind stress term in force balance (x) 
    306297      IF ( iom_use( "vtau_ice" ) ) CALL iom_put( "vtau_ice"     ,  vtau_ice(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Wind stress term in force balance (y) 
    307  
    308  
    309       IF ( iom_use( "icestr"   ) ) CALL iom_put( "icestr"      ,   strength(:,:)            * zswi(:,:) + zmiss(:,:) )   ! Ice strength 
    310298 
    311299      !-------------------------------- 
Note: See TracChangeset for help on using the changeset viewer.