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 4332 for branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90 – NEMO

Ignore:
Timestamp:
2013-12-11T15:38:42+01:00 (10 years ago)
Author:
clem
Message:

update LIM3 to fix remaining bugs. Now working in global and regional config.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2013/dev_r4028_CNRS_LIM3/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90

    r4072 r4332  
    4343   PUBLIC   lim_itd_me_alloc        ! called by iceini.F90 
    4444 
    45    REAL(wp) ::   epsi11 = 1.e-11_wp   ! constant values 
     45   REAL(wp) ::   epsi20 = 1.e-20_wp   ! constant values 
    4646   REAL(wp) ::   epsi10 = 1.e-10_wp   ! constant values 
    4747   REAL(wp) ::   epsi06 = 1.e-06_wp   ! constant values 
     
    5050   ! Variables shared among ridging subroutines 
    5151   !----------------------------------------------------------------------- 
    52    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   asum     ! sum of total ice and open water area 
    53    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   aksum    ! ratio of area removed to area ridged 
    54    ! 
    55    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   athorn   ! participation function; fraction of ridging/ 
    56    !                                                           !  closing associated w/ category n 
    57    ! 
    58    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hrmin    ! minimum ridge thickness 
    59    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hrmax    ! maximum ridge thickness 
    60    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   hraft    ! thickness of rafted ice 
    61    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   krdg     ! mean ridge thickness/thickness of ridging ice  
    62    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   aridge   ! participating ice ridging 
    63    REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   araft    ! participating ice rafting 
     52   REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   asum     ! sum of total ice and open water area 
     53   REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   aksum    ! ratio of area removed to area ridged 
     54   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   athorn   ! participation function; fraction of ridging/ 
     55   !                                                     ! closing associated w/ category n 
     56   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   hrmin    ! minimum ridge thickness 
     57   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   hrmax    ! maximum ridge thickness 
     58   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   hraft    ! thickness of rafted ice 
     59   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   krdg     ! mean ridge thickness/thickness of ridging ice  
     60   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   aridge   ! participating ice ridging 
     61   REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   araft    ! participating ice rafting 
    6462 
    6563   REAL(wp), PARAMETER ::   krdgmin = 1.1_wp    ! min ridge thickness multiplier 
    6664   REAL(wp), PARAMETER ::   kraft   = 2.0_wp    ! rafting multipliyer 
    67    REAL(wp), PARAMETER ::   kamax   = 1.0 
     65   REAL(wp), PARAMETER ::   kamax   = 1.0_wp    ! max of ice area authorized (clem: scheme is not stable if kamax <= 0.99) 
    6866 
    6967   REAL(wp) ::   Cp                             !  
     
    183181      ! 1) Thickness categories boundaries, ice / o.w. concentrations, init_ons 
    184182      !-----------------------------------------------------------------------------! 
    185       ! Set hi_max(ncat) to a big value to ensure that all ridged ice is thinner than hi_max(ncat). 
    186  
    187       hi_max(jpl) = 999.99 
    188  
    189183      Cp = 0.5 * grav * (rau0-rhoic) * rhoic / rau0                ! proport const for PE 
    190184      ! 
     
    265259               ! Reduce the closing rate if more than 100% of the open water  
    266260               ! would be removed.  Reduce the opening rate proportionately. 
    267                IF ( ato_i(ji,jj) .GT. epsi11 .AND. athorn(ji,jj,0) .GT. 0.0 ) THEN 
     261               IF ( ato_i(ji,jj) .GT. epsi10 .AND. athorn(ji,jj,0) .GT. 0.0 ) THEN 
    268262                  w1 = athorn(ji,jj,0) * closing_gross(ji,jj) * rdt_ice 
    269263                  IF ( w1 .GT. ato_i(ji,jj)) THEN 
     
    285279            DO jj = 1, jpj 
    286280               DO ji = 1, jpi 
    287                   IF ( a_i(ji,jj,jl) > epsi11 .AND. athorn(ji,jj,jl) > 0._wp )THEN 
     281                  IF ( a_i(ji,jj,jl) > epsi10 .AND. athorn(ji,jj,jl) > 0._wp )THEN 
    288282                     w1 = athorn(ji,jj,jl) * closing_gross(ji,jj) * rdt_ice 
    289283                     IF ( w1  >  a_i(ji,jj,jl) ) THEN 
     
    316310         DO jj = 1, jpj 
    317311            DO ji = 1, jpi 
    318                IF (ABS(asum(ji,jj) - kamax ) .LT. epsi11) THEN 
     312               IF (ABS(asum(ji,jj) - kamax ) .LT. epsi10) THEN 
    319313                  closing_net(ji,jj) = 0._wp 
    320314                  opning     (ji,jj) = 0._wp 
     
    358352         DO ji = 1, jpi 
    359353 
    360             IF(ABS(asum(ji,jj) - kamax) > epsi11 ) asum_error = .true. 
     354            IF(ABS(asum(ji,jj) - kamax) > epsi10 ) asum_error = .true. 
    361355 
    362356            dardg1dt(ji,jj) = dardg1dt(ji,jj) * r1_rdtice 
     
    377371      DO jj = 1, jpj 
    378372         DO ji = 1, jpi 
    379             IF( ABS( asum(ji,jj) - kamax)  >  epsi11 ) THEN   ! there is a bug 
     373            IF( ABS( asum(ji,jj) - kamax)  >  epsi10 ) THEN   ! there is a bug 
    380374               WRITE(numout,*) ' ' 
    381375               WRITE(numout,*) ' ALERTE : Ridging error: total area = ', asum(ji,jj) 
     
    400394 
    401395      !-----------------------------------------------------------------------------! 
    402       ! 6) Updating state variables and trend terms 
     396      ! 6) Updating state variables and trend terms (done in limupdate) 
    403397      !-----------------------------------------------------------------------------! 
    404  
    405398      CALL lim_var_glo2eqv 
    406399      CALL lim_itd_me_zapsmall 
     
    419412         END DO 
    420413      END DO 
    421  
    422       !----------------- 
    423       !  Trend terms 
    424       !----------------- 
    425  
    426       d_u_ice_dyn(:,:)     = u_ice(:,:)     - old_u_ice(:,:) 
    427       d_v_ice_dyn(:,:)     = v_ice(:,:)     - old_v_ice(:,:) 
    428       d_a_i_trp  (:,:,:)   = a_i  (:,:,:)   - old_a_i  (:,:,:) 
    429       d_v_s_trp  (:,:,:)   = v_s  (:,:,:)   - old_v_s  (:,:,:)   
    430       d_v_i_trp  (:,:,:)   = v_i  (:,:,:)   - old_v_i  (:,:,:)    
    431       d_e_s_trp  (:,:,:,:) = e_s  (:,:,:,:) - old_e_s  (:,:,:,:)   
    432       d_e_i_trp  (:,:,:,:) = e_i  (:,:,:,:) - old_e_i  (:,:,:,:) 
    433       d_oa_i_trp (:,:,:)   = oa_i (:,:,:)   - old_oa_i (:,:,:) 
    434       d_smv_i_trp(:,:,:)   = 0._wp 
    435       IF(  num_sal == 2  )   d_smv_i_trp(:,:,:)  = smv_i(:,:,:) - old_smv_i(:,:,:) 
    436414 
    437415      IF(ln_ctl) THEN     ! Control print 
     
    491469      ! ------------------------------- 
    492470 
    493       !-------------------------! 
    494       ! Back to initial values 
    495       !-------------------------! 
    496  
    497       ! update of fields will be made later in lim update 
    498       u_ice(:,:)   = old_u_ice(:,:) 
    499       v_ice(:,:)   = old_v_ice(:,:) 
    500       a_i(:,:,:)   = old_a_i(:,:,:) 
    501       v_s(:,:,:)   = old_v_s(:,:,:) 
    502       v_i(:,:,:)   = old_v_i(:,:,:) 
    503       e_s(:,:,:,:) = old_e_s(:,:,:,:) 
    504       e_i(:,:,:,:) = old_e_i(:,:,:,:) 
    505       oa_i(:,:,:)  = old_oa_i(:,:,:) 
    506       IF(  num_sal == 2  )   smv_i(:,:,:) = old_smv_i(:,:,:) 
    507  
    508       !----------------------------------------------------! 
    509       ! Advection of ice in a free cell, newly ridged ice 
    510       !----------------------------------------------------! 
    511  
    512       ! to allow for thermodynamics to melt new ice 
    513       ! we immediately advect ice in free cells 
    514  
    515       ! heat content has to be corrected before ice volume 
    516 !clem@order 
    517 !      DO jl = 1, jpl 
    518 !         DO jk = 1, nlay_i 
    519 !            DO jj = 1, jpj 
    520 !               DO ji = 1, jpi 
    521 !                  IF ( ( old_v_i(ji,jj,jl) < epsi06 ) .AND. & 
    522 !                     ( d_v_i_trp(ji,jj,jl) > epsi06 ) ) THEN 
    523 !                     old_e_i(ji,jj,jk,jl)   = d_e_i_trp(ji,jj,jk,jl) 
    524 !                     d_e_i_trp(ji,jj,jk,jl) = 0._wp 
    525 !                  ENDIF 
    526 !               END DO 
    527 !            END DO 
    528 !         END DO 
    529 !      END DO 
    530 ! 
    531 !      DO jl = 1, jpl 
    532 !         DO jj = 1, jpj 
    533 !            DO ji = 1, jpi 
    534 !               IF(  old_v_i  (ji,jj,jl) < epsi06  .AND. & 
    535 !                    d_v_i_trp(ji,jj,jl) > epsi06    ) THEN 
    536 !                  old_v_i   (ji,jj,jl)   = d_v_i_trp(ji,jj,jl) 
    537 !                  d_v_i_trp (ji,jj,jl)   = 0._wp 
    538 !                  old_a_i   (ji,jj,jl)   = d_a_i_trp(ji,jj,jl) 
    539 !                  d_a_i_trp (ji,jj,jl)   = 0._wp 
    540 !                  old_v_s   (ji,jj,jl)   = d_v_s_trp(ji,jj,jl) 
    541 !                  d_v_s_trp (ji,jj,jl)   = 0._wp 
    542 !                  old_e_s   (ji,jj,1,jl) = d_e_s_trp(ji,jj,1,jl) 
    543 !                  d_e_s_trp (ji,jj,1,jl) = 0._wp 
    544 !                  old_oa_i  (ji,jj,jl)   = d_oa_i_trp(ji,jj,jl) 
    545 !                  d_oa_i_trp(ji,jj,jl)   = 0._wp 
    546 !                  IF(  num_sal == 2  )   old_smv_i(ji,jj,jl) = d_smv_i_trp(ji,jj,jl) 
    547 !                  d_smv_i_trp(ji,jj,jl)  = 0._wp 
    548 !               ENDIF 
    549 !            END DO 
    550 !         END DO 
    551 !      END DO 
    552 !clem@order 
    553471      ENDIF  ! ln_limdyn=.true. 
    554472      ! 
     
    605523               DO ji = 1, jpi 
    606524                  ! 
    607                   IF(  a_i(ji,jj,jl)    > epsi11  .AND.     & 
    608                        athorn(ji,jj,jl) > 0._wp   ) THEN 
     525                  IF( a_i(ji,jj,jl) > epsi10 .AND. athorn(ji,jj,jl) > 0._wp ) THEN 
    609526                     hi = v_i(ji,jj,jl) / a_i(ji,jj,jl) 
    610527                     !---------------------------- 
     
    624541                        * z1_3 * (hrmax(ji,jj,jl)**3 - hrmin(ji,jj,jl)**3) / ( hrmax(ji,jj,jl)-hrmin(ji,jj,jl) )    
    625542!!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...                     
    626                   ENDIF            ! aicen > epsi11 
     543                  ENDIF            ! aicen > epsi10 
    627544                  ! 
    628545               END DO ! ji 
     
    681598         DO jj = 2, jpj - 1 
    682599            DO ji = 2, jpi - 1 
    683                IF ( ( asum(ji,jj) - ato_i(ji,jj) ) .GT. epsi11) THEN ! ice is 
     600               IF ( ( asum(ji,jj) - ato_i(ji,jj) ) .GT. epsi10) THEN ! ice is 
    684601                  ! present 
    685602                  zworka(ji,jj) = 4.0 * strength(ji,jj)              & 
     
    721638         DO jj = 1, jpj - 1 
    722639            DO ji = 1, jpi - 1 
    723                IF ( ( asum(ji,jj) - ato_i(ji,jj) ) .GT. epsi11) THEN       ! ice is present 
     640               IF ( ( asum(ji,jj) - ato_i(ji,jj) ) .GT. epsi10) THEN       ! ice is present 
    724641                  numts_rm = 1 ! number of time steps for the running mean 
    725642                  IF ( strp1(ji,jj) .GT. 0.0 ) numts_rm = numts_rm + 1 
     
    794711      DO jj = 1, jpj 
    795712         DO ji = 1, jpi 
    796             IF( ato_i(ji,jj) > epsi11 ) THEN   ;   Gsum(ji,jj,0) = ato_i(ji,jj) 
     713            IF( ato_i(ji,jj) > epsi10 ) THEN   ;   Gsum(ji,jj,0) = ato_i(ji,jj) 
    797714            ELSE                               ;   Gsum(ji,jj,0) = 0._wp 
    798715            ENDIF 
     
    804721         DO jj = 1, jpj  
    805722            DO ji = 1, jpi 
    806                IF( a_i(ji,jj,jl) .GT. epsi11 ) THEN   ;   Gsum(ji,jj,jl) = Gsum(ji,jj,jl-1) + a_i(ji,jj,jl) 
     723               IF( a_i(ji,jj,jl) .GT. epsi10 ) THEN   ;   Gsum(ji,jj,jl) = Gsum(ji,jj,jl-1) + a_i(ji,jj,jl) 
    807724               ELSE                                   ;   Gsum(ji,jj,jl) = Gsum(ji,jj,jl-1) 
    808725               ENDIF 
     
    887804      IF ( raftswi == 1 ) THEN 
    888805 
    889          IF( MAXVAL(aridge + araft - athorn(:,:,1:jpl)) .GT. epsi11 ) THEN 
     806         IF( MAXVAL(aridge + araft - athorn(:,:,1:jpl)) .GT. epsi10 ) THEN 
    890807            DO jl = 1, jpl 
    891808               DO jj = 1, jpj 
    892809                  DO ji = 1, jpi 
    893                      IF ( aridge(ji,jj,jl) + araft(ji,jj,jl) - athorn(ji,jj,jl) .GT. & 
    894                         epsi11 ) THEN 
     810                     IF ( aridge(ji,jj,jl) + araft(ji,jj,jl) - athorn(ji,jj,jl) .GT. epsi10 ) THEN 
    895811                        WRITE(numout,*) ' ALERTE 96 : wrong participation function ... ' 
    896812                        WRITE(numout,*) ' ji, jj, jl : ', ji, jj, jl 
     
    938854            DO ji = 1, jpi 
    939855 
    940                IF (a_i(ji,jj,jl) .GT. epsi11 .AND. athorn(ji,jj,jl) .GT. 0.0 ) THEN 
     856               IF (a_i(ji,jj,jl) .GT. epsi10 .AND. athorn(ji,jj,jl) .GT. 0.0 ) THEN 
    941857                  hi = v_i(ji,jj,jl) / a_i(ji,jj,jl) 
    942858                  hrmean          = MAX(SQRT(Hstar*hi), hi*krdgmin) 
     
    992908      INTEGER ::   ij                ! horizontal index, combines i and j loops 
    993909      INTEGER ::   icells            ! number of cells with aicen > puny 
    994       REAL(wp) ::   zeps, zindb, zsrdg2   ! local scalar 
     910      REAL(wp) ::   zindb, zsrdg2   ! local scalar 
    995911      REAL(wp) ::   hL, hR, farea, zdummy, zdummy0, ztmelts    ! left and right limits of integration 
    996912 
     
    1044960 
    1045961      IF( con_i ) THEN 
    1046          CALL lim_column_sum (jpl,   v_i, vice_init ) 
    1047          WRITE(numout,*) ' vice_init  : ', vice_init(jiindx,jjindx) 
     962         CALL lim_column_sum        (jpl,    v_i,       vice_init ) 
    1048963         CALL lim_column_sum_energy (jpl, nlay_i,  e_i, eice_init ) 
    1049          WRITE(numout,*) ' eice_init  : ', eice_init(jiindx,jjindx) 
     964         DO ji = mi0(jiindx), mi1(jiindx) 
     965            DO jj = mj0(jjindx), mj1(jjindx) 
     966               WRITE(numout,*) ' vice_init  : ', vice_init(ji,jj) 
     967               WRITE(numout,*) ' eice_init  : ', eice_init(ji,jj) 
     968            END DO 
     969         END DO 
    1050970      ENDIF 
    1051  
    1052       zeps   = 1.e-20_wp 
    1053971 
    1054972      !------------------------------------------------------------------------------- 
     
    1062980            ato_i(ji,jj) = ato_i(ji,jj) - athorn(ji,jj,0) * closing_gross(ji,jj) * rdt_ice        & 
    1063981               &                        + opning(ji,jj)                          * rdt_ice 
    1064             IF( ato_i(ji,jj) < -epsi11 ) THEN 
     982            IF( ato_i(ji,jj) < -epsi10 ) THEN 
    1065983               neg_ato_i = .TRUE. 
    1066984            ELSEIF( ato_i(ji,jj) < 0._wp ) THEN    ! roundoff error 
     
    1074992         DO jj = 1, jpj  
    1075993            DO ji = 1, jpi 
    1076                IF( ato_i(ji,jj) < -epsi11 ) THEN  
     994               IF( ato_i(ji,jj) < -epsi10 ) THEN  
    1077995                  WRITE(numout,*) ''   
    1078996                  WRITE(numout,*) 'Ridging error: ato_i < 0' 
    1079997                  WRITE(numout,*) 'ato_i : ', ato_i(ji,jj) 
    1080                ENDIF               ! ato_i < -epsi11 
     998               ENDIF               ! ato_i < -epsi10 
    1081999            END DO 
    10821000         END DO 
     
    11201038         DO jj = 1, jpj 
    11211039            DO ji = 1, jpi 
    1122                IF( aicen_init(ji,jj,jl1)  >  epsi11  .AND.  athorn(ji,jj,jl1)  >  0._wp       & 
    1123                   .AND. closing_gross(ji,jj) > 0._wp ) THEN 
     1040               IF( aicen_init(ji,jj,jl1) > epsi10 .AND. athorn(ji,jj,jl1) > 0._wp  & 
     1041                  &   .AND. closing_gross(ji,jj) > 0._wp ) THEN 
    11241042                  icells = icells + 1 
    11251043                  indxi(icells) = ji 
     
    11581076            afrft(ji,jj) = arft1(ji,jj) / aicen_init(ji,jj,jl1) !rafting 
    11591077 
    1160             IF (afrac(ji,jj) > kamax + epsi11) THEN  !riging 
     1078            IF (afrac(ji,jj) > kamax + epsi10) THEN  !riging 
    11611079               large_afrac = .true. 
    11621080            ELSEIF (afrac(ji,jj) > kamax) THEN  ! roundoff error 
    11631081               afrac(ji,jj) = kamax 
    11641082            ENDIF 
    1165             IF (afrft(ji,jj) > kamax + epsi11) THEN !rafting 
     1083            IF (afrft(ji,jj) > kamax + epsi10) THEN !rafting 
    11661084               large_afrft = .true. 
    11671085            ELSEIF (afrft(ji,jj) > kamax) THEN  ! roundoff error 
     
    12221140            dardg1dt   (ji,jj) = dardg1dt(ji,jj) + ardg1(ji,jj) + arft1(ji,jj) 
    12231141            dardg2dt   (ji,jj) = dardg2dt(ji,jj) + ardg2(ji,jj) + arft2(ji,jj) 
    1224             !clem diag_dyn_gr(ji,jj) = diag_dyn_gr(ji,jj) + ( vrdg2(ji,jj) + virft(ji,jj) ) * r1_rdtice 
    12251142            opening    (ji,jj) = opening (ji,jj) + opning(ji,jj) * rdt_ice 
    12261143 
     
    12761193 
    12771194               ! corrected sea water salinity 
    1278                zindb  = MAX( 0._wp , SIGN( 1._wp , vsw(ji,jj) - zeps ) ) 
    1279                zdummy = zindb * ( srdg1(ji,jj) - srdg2(ji,jj) ) / MAX( ridge_por * vsw(ji,jj), zeps ) 
     1195               zindb  = MAX( 0._wp , SIGN( 1._wp , vsw(ji,jj) - epsi20 ) ) 
     1196               zdummy = zindb * ( srdg1(ji,jj) - srdg2(ji,jj) ) / MAX( ridge_por * vsw(ji,jj), epsi20 ) 
    12801197 
    12811198               ztmelts          = - tmut * zdummy + rtt 
     
    13121229               ji = indxi(ij) 
    13131230               jj = indxj(ij) 
    1314                IF( afrac(ji,jj) > kamax + epsi11 ) THEN  
     1231               IF( afrac(ji,jj) > kamax + epsi10 ) THEN  
    13151232                  WRITE(numout,*) '' 
    13161233                  WRITE(numout,*) ' ardg > a_i' 
     
    13241241               ji = indxi(ij) 
    13251242               jj = indxj(ij) 
    1326                IF( afrft(ji,jj) > kamax + epsi11 ) THEN  
     1243               IF( afrft(ji,jj) > kamax + epsi10 ) THEN  
    13271244                  WRITE(numout,*) '' 
    13281245                  WRITE(numout,*) ' arft > a_i' 
     
    14241341         fieldid = ' v_i : limitd_me ' 
    14251342         CALL lim_cons_check (vice_init, vice_final, 1.0e-6, fieldid)  
    1426          WRITE(numout,*) ' vice_init  : ', vice_init(jiindx,jjindx) 
    1427          WRITE(numout,*) ' vice_final : ', vice_final(jiindx,jjindx) 
    14281343 
    14291344         CALL lim_column_sum_energy (jpl, nlay_i,  e_i, eice_final ) 
    14301345         fieldid = ' e_i : limitd_me ' 
    14311346         CALL lim_cons_check (eice_init, eice_final, 1.0e-2, fieldid)  
    1432          WRITE(numout,*) ' eice_init  : ', eice_init(jiindx,jjindx) 
    1433          WRITE(numout,*) ' eice_final : ', eice_final(jiindx,jjindx) 
     1347 
     1348         DO ji = mi0(jiindx), mi1(jiindx) 
     1349            DO jj = mj0(jjindx), mj1(jjindx) 
     1350               WRITE(numout,*) ' vice_init  : ', vice_init (ji,jj) 
     1351               WRITE(numout,*) ' vice_final : ', vice_final(ji,jj) 
     1352               WRITE(numout,*) ' eice_init  : ', eice_init (ji,jj) 
     1353               WRITE(numout,*) ' eice_final : ', eice_final(ji,jj) 
     1354            END DO 
     1355         END DO 
    14341356      ENDIF 
    14351357      ! 
     
    15331455 
    15341456      REAL(wp), POINTER, DIMENSION(:,:) ::   zmask   ! 2D workspace 
    1535        
     1457      REAL(wp)                          ::   zmask_glo 
    15361458!!gm      REAL(wp) ::   xtmp      ! temporary variable 
    15371459      !!------------------------------------------------------------------- 
     
    15451467         ! Abort model in case of negative area. 
    15461468         !----------------------------------------------------------------- 
    1547          IF( MINVAL(a_i(:,:,jl)) .LT. -epsi11 .AND. ln_nicep ) THEN 
    1548             DO jj = 1, jpj 
    1549                DO ji = 1, jpi 
    1550                   IF ( a_i(ji,jj,jl) .LT. -epsi11 ) THEN 
    1551                      WRITE (numout,*) ' ALERTE 98 '  
    1552                      WRITE (numout,*) ' Negative ice area: ji, jj, jl: ', ji, jj,jl 
    1553                      WRITE (numout,*) ' a_i    ', a_i(ji,jj,jl) 
    1554                   ENDIF 
    1555                END DO 
    1556             END DO 
    1557          ENDIF 
    1558  
    15591469         icells = 0 
    1560          zmask  = 0._wp 
     1470         zmask(:,:)  = 0._wp 
    15611471         DO jj = 1, jpj 
    15621472            DO ji = 1, jpi 
    1563                IF(  ( a_i(ji,jj,jl) .GE. -epsi11 .AND. a_i(ji,jj,jl) .LT. 0._wp   ) .OR.   & 
    1564                   & ( a_i(ji,jj,jl) .GT. 0._wp   .AND. a_i(ji,jj,jl) .LE. 1.0e-11 ) .OR.   & 
    1565                   & ( v_i(ji,jj,jl)  ==  0._wp   .AND. a_i(ji,jj,jl) .GT. 0._wp   ) .OR.   & 
    1566                   & ( v_i(ji,jj,jl) .GT. 0._wp   .AND. v_i(ji,jj,jl) .LT. 1.e-12  )      )   zmask(ji,jj) = 1._wp 
    1567             END DO 
    1568          END DO 
    1569          IF( ln_nicep )   WRITE(numout,*) SUM(zmask), ' cells of ice zapped in the ocean ' 
     1473               IF(  ( a_i(ji,jj,jl) >= -epsi10 .AND. a_i(ji,jj,jl) <  0._wp   ) .OR.   & 
     1474                  & ( a_i(ji,jj,jl) >  0._wp   .AND. a_i(ji,jj,jl) <= epsi10  ) .OR.   & 
     1475                  & ( v_i(ji,jj,jl) == 0._wp   .AND. a_i(ji,jj,jl) >  0._wp   ) .OR.   & 
     1476                  & ( v_i(ji,jj,jl) >  0._wp   .AND. v_i(ji,jj,jl) <= epsi10  ) )   zmask(ji,jj) = 1._wp 
     1477            END DO 
     1478         END DO 
     1479         zmask_glo = glob_sum(zmask) 
     1480         !IF( ln_nicep .AND. lwp ) WRITE(numout,*) zmask_glo, ' cells of ice zapped in the ocean ' 
    15701481 
    15711482         !----------------------------------------------------------------- 
     
    15791490!!gm                  xtmp = xtmp * unit_fac 
    15801491                  ! fheat_res(ji,jj) = fheat_res(ji,jj) - xtmp 
    1581                   e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * ( 1 - zmask(ji,jj) ) 
     1492                  e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * ( 1._wp - zmask(ji,jj) ) 
    15821493               END DO 
    15831494            END DO 
     
    16011512!!gm               xtmp = ( rhosn*cpic*( rtt-t_s(ji,jj,1,jl) ) + rhosn*lfus ) * r1_rdtice !RB   ??????? 
    16021513 
    1603                t_s(ji,jj,1,jl) = rtt * zmask(ji,jj) + t_s(ji,jj,1,jl) * ( 1 - zmask(ji,jj) ) 
     1514               t_s(ji,jj,1,jl) = rtt * zmask(ji,jj) + t_s(ji,jj,1,jl) * ( 1._wp - zmask(ji,jj) ) 
    16041515 
    16051516               !----------------------------------------------------------------- 
     
    16151526 
    16161527               ato_i(ji,jj)    = a_i  (ji,jj,jl) *       zmask(ji,jj)   + ato_i(ji,jj) 
    1617                a_i  (ji,jj,jl) = a_i  (ji,jj,jl) * ( 1 - zmask(ji,jj) ) 
    1618                v_i  (ji,jj,jl) = v_i  (ji,jj,jl) * ( 1 - zmask(ji,jj) ) 
    1619                v_s  (ji,jj,jl) = v_s  (ji,jj,jl) * ( 1 - zmask(ji,jj) ) 
    1620                t_su (ji,jj,jl) = t_su (ji,jj,jl) * ( 1 - zmask(ji,jj) ) + t_bo(ji,jj) * zmask(ji,jj) 
    1621                oa_i (ji,jj,jl) = oa_i (ji,jj,jl) * ( 1 - zmask(ji,jj) ) 
    1622                smv_i(ji,jj,jl) = smv_i(ji,jj,jl) * ( 1 - zmask(ji,jj) ) 
     1528               a_i  (ji,jj,jl) = a_i  (ji,jj,jl) * ( 1._wp - zmask(ji,jj) ) 
     1529               v_i  (ji,jj,jl) = v_i  (ji,jj,jl) * ( 1._wp - zmask(ji,jj) ) 
     1530               v_s  (ji,jj,jl) = v_s  (ji,jj,jl) * ( 1._wp - zmask(ji,jj) ) 
     1531               t_su (ji,jj,jl) = t_su (ji,jj,jl) * ( 1._wp - zmask(ji,jj) ) + t_bo(ji,jj) * zmask(ji,jj) 
     1532               oa_i (ji,jj,jl) = oa_i (ji,jj,jl) * ( 1._wp - zmask(ji,jj) ) 
     1533               smv_i(ji,jj,jl) = smv_i(ji,jj,jl) * ( 1._wp - zmask(ji,jj) ) 
    16231534               ! 
    16241535            END DO 
Note: See TracChangeset for help on using the changeset viewer.