Changeset 5080


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

LIM3 cleaning again (almost done)

Location:
branches/2015/dev_r5044_CNRS_LIM3CLEAN/NEMOGCM/NEMO
Files:
6 edited

Legend:

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

    r5078 r5080  
    424424 
    425425      ii = ii + 1 
    426       ALLOCATE( sist   (jpi,jpj) , icethi (jpi,jpj) , t_bo   (jpi,jpj) ,      & 
    427          &      frld   (jpi,jpj) , pfrld  (jpi,jpj) , phicif (jpi,jpj) ,      & 
    428          &      wfx_snw(jpi,jpj) , wfx_ice(jpi,jpj) , wfx_sub(jpi,jpj) ,    & 
     426      ALLOCATE( sist   (jpi,jpj) , icethi (jpi,jpj) , t_bo   (jpi,jpj) ,                        & 
     427         &      frld   (jpi,jpj) , pfrld  (jpi,jpj) , phicif (jpi,jpj) ,                        & 
     428         &      wfx_snw(jpi,jpj) , wfx_ice(jpi,jpj) , wfx_sub(jpi,jpj) ,                        & 
    429429         &      wfx_bog(jpi,jpj) , wfx_dyn(jpi,jpj) , wfx_bom(jpi,jpj) , wfx_sum(jpi,jpj) ,     & 
    430          &      wfx_res(jpi,jpj) , wfx_sni(jpi,jpj) , wfx_opw(jpi,jpj) , wfx_spr(jpi,jpj) ,  qlead  (jpi,jpj) ,     & 
    431          &      afx_tot(jpi,jpj) , afx_thd(jpi,jpj),  afx_dyn(jpi,jpj) , & 
    432          &      fhtur  (jpi,jpj) , ftr_ice(jpi,jpj,jpl),    & 
    433          &      sfx_res(jpi,jpj) , sfx_bri(jpi,jpj) , sfx_dyn(jpi,jpj) ,      & 
    434          &      sfx_bog(jpi,jpj) , sfx_bom(jpi,jpj) , sfx_sum(jpi,jpj) , sfx_sni(jpi,jpj) , sfx_opw(jpi,jpj) ,   & 
     430         &      wfx_res(jpi,jpj) , wfx_sni(jpi,jpj) , wfx_opw(jpi,jpj) , wfx_spr(jpi,jpj) ,     & 
     431         &      afx_tot(jpi,jpj) , afx_thd(jpi,jpj),  afx_dyn(jpi,jpj) ,                        & 
     432         &      fhtur  (jpi,jpj) , ftr_ice(jpi,jpj,jpl), qlead  (jpi,jpj) ,                     & 
     433         &      sfx_res(jpi,jpj) , sfx_bri(jpi,jpj) , sfx_dyn(jpi,jpj) ,                        & 
     434         &      sfx_bog(jpi,jpj) , sfx_bom(jpi,jpj) , sfx_sum(jpi,jpj) , sfx_sni(jpi,jpj) , sfx_opw(jpi,jpj) ,    & 
    435435         &      hfx_res(jpi,jpj) , hfx_snw(jpi,jpj) , hfx_sub(jpi,jpj) , hfx_err(jpi,jpj) , hfx_err_rem(jpi,jpj), & 
    436          &      hfx_in (jpi,jpj) , hfx_out(jpi,jpj) , fhld(jpi,jpj) ,  & 
    437          &      hfx_sum(jpi,jpj) , hfx_bom(jpi,jpj) , hfx_bog(jpi,jpj) , hfx_dif(jpi,jpj) , hfx_opw(jpi,jpj) , & 
     436         &      hfx_in (jpi,jpj) , hfx_out(jpi,jpj) , fhld(jpi,jpj) ,                           & 
     437         &      hfx_sum(jpi,jpj) , hfx_bom(jpi,jpj) , hfx_bog(jpi,jpj) , hfx_dif(jpi,jpj) , hfx_opw(jpi,jpj) ,    & 
    438438         &      hfx_thd(jpi,jpj) , hfx_dyn(jpi,jpj) , hfx_spr(jpi,jpj) ,  STAT=ierr(ii) ) 
    439439 
     
    450450         &      bv_i (jpi,jpj) , smt_i(jpi,jpj)                                   , STAT=ierr(ii) ) 
    451451      ii = ii + 1 
    452       ALLOCATE( t_s(jpi,jpj,nlay_s,jpl) ,                            & 
    453          &      e_s(jpi,jpj,nlay_s,jpl) , STAT=ierr(ii) ) 
     452      ALLOCATE( t_s(jpi,jpj,nlay_s,jpl) , e_s(jpi,jpj,nlay_s,jpl) , STAT=ierr(ii) ) 
    454453      ii = ii + 1 
    455454      ALLOCATE( t_i(jpi,jpj,nlay_i+1,jpl) , e_i(jpi,jpj,nlay_i+1,jpl) , s_i(jpi,jpj,nlay_i+1,jpl) , STAT=ierr(ii) ) 
     
    475474      ii = ii + 1 
    476475      ALLOCATE( v_s_b  (jpi,jpj,jpl) , v_i_b  (jpi,jpj,jpl) , e_s_b(jpi,jpj,nlay_s,jpl) ,     & 
    477          &      a_i_b  (jpi,jpj,jpl) , smv_i_b(jpi,jpj,jpl) , e_i_b(jpi,jpj,nlay_i+1 ,jpl) ,     & 
    478          &      oa_i_b (jpi,jpj,jpl)                                                        ,     & 
    479          &      u_ice_b(jpi,jpj)     , v_ice_b(jpi,jpj)                                   , STAT=ierr(ii) ) 
     476         &      a_i_b  (jpi,jpj,jpl) , smv_i_b(jpi,jpj,jpl) , e_i_b(jpi,jpj,nlay_i+1 ,jpl) ,  & 
     477         &      oa_i_b (jpi,jpj,jpl) , u_ice_b(jpi,jpj)     , v_ice_b(jpi,jpj)             , STAT=ierr(ii) ) 
    480478       
    481479      ! * Ice thickness distribution variables 
  • branches/2015/dev_r5044_CNRS_LIM3CLEAN/NEMOGCM/NEMO/LIM_SRC_3/limadv.F90

    r5067 r5080  
    6363      !!  
    6464      INTEGER  ::   ji, jj                               ! dummy loop indices 
    65       REAL(wp) ::   zs1max, zrdt, zslpmax, ztemp, zin0   ! local scalars 
     65      REAL(wp) ::   zs1max, zrdt, zslpmax, ztemp         ! local scalars 
    6666      REAL(wp) ::   zs1new, zalf , zalfq , zbt           !   -      - 
    6767      REAL(wp) ::   zs2new, zalf1, zalf1q, zbt1          !   -      - 
     
    8585            zs2new  = MIN(  2.0 * zslpmax - 0.3334 * ABS( zs1new ),      & 
    8686               &            MAX( ABS( zs1new ) - zslpmax, psxx(ji,jj) )  ) 
    87             zin0    = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask 
     87            rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask 
    8888 
    8989            ps0 (ji,jj) = zslpmax   
    90             psx (ji,jj) = zs1new      * zin0 
    91             psxx(ji,jj) = zs2new      * zin0 
    92             psy (ji,jj) = psy (ji,jj) * zin0 
    93             psyy(ji,jj) = psyy(ji,jj) * zin0 
    94             psxy(ji,jj) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj) ) ) * zin0 
     90            psx (ji,jj) = zs1new      * rswitch 
     91            psxx(ji,jj) = zs2new      * rswitch 
     92            psy (ji,jj) = psy (ji,jj) * rswitch 
     93            psyy(ji,jj) = psyy(ji,jj) * rswitch 
     94            psxy(ji,jj) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj) ) ) * rswitch 
    9595         END DO 
    9696      END DO 
     
    248248      !! 
    249249      INTEGER  ::   ji, jj                               ! dummy loop indices 
    250       REAL(wp) ::   zs1max, zrdt, zslpmax, ztemp, zin0   ! temporary scalars 
     250      REAL(wp) ::   zs1max, zrdt, zslpmax, ztemp         ! temporary scalars 
    251251      REAL(wp) ::   zs1new, zalf , zalfq , zbt           !    -         - 
    252252      REAL(wp) ::   zs2new, zalf1, zalf1q, zbt1          !    -         - 
     
    270270            zs2new  = MIN(  ( 2.0 * zslpmax - 0.3334 * ABS( zs1new ) ),   & 
    271271               &             MAX( ABS( zs1new )-zslpmax, psyy(ji,jj) )  ) 
    272             zin0    = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask 
     272            rswitch = ( 1.0 - MAX( 0._wp, SIGN( 1._wp, -zslpmax) ) ) * tmask(ji,jj,1)   ! Case of empty boxes & Apply mask 
    273273            ! 
    274274            ps0 (ji,jj) = zslpmax   
    275             psx (ji,jj) = psx (ji,jj) * zin0 
    276             psxx(ji,jj) = psxx(ji,jj) * zin0 
    277             psy (ji,jj) = zs1new * zin0 
    278             psyy(ji,jj) = zs2new * zin0 
    279             psxy(ji,jj) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj) ) ) * zin0 
     275            psx (ji,jj) = psx (ji,jj) * rswitch 
     276            psxx(ji,jj) = psxx(ji,jj) * rswitch 
     277            psy (ji,jj) = zs1new * rswitch 
     278            psyy(ji,jj) = zs2new * rswitch 
     279            psxy(ji,jj) = MIN( zslpmax, MAX( -zslpmax, psxy(ji,jj) ) ) * rswitch 
    280280         END DO 
    281281      END DO 
  • branches/2015/dev_r5044_CNRS_LIM3CLEAN/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90

    r5070 r5080  
    103103            DO jj = 1, jpj 
    104104               zswitch(jj) = SUM( 1.0 - at_i(:,jj) )   ! = REAL(jpj) if ocean everywhere on a j-line 
    105                zmsk(jj) = SUM( tmask(:,jj,1)    )   ! = 0         if land  everywhere on a j-line 
     105               zmsk   (jj) = SUM( tmask(:,jj,1)    )   ! = 0         if land  everywhere on a j-line 
    106106            END DO 
    107107 
  • branches/2015/dev_r5044_CNRS_LIM3CLEAN/NEMOGCM/NEMO/LIM_SRC_3/limhdf.F90

    r5070 r5080  
    2626   PRIVATE 
    2727 
    28    PUBLIC   lim_hdf     ! called by lim_tra 
     28   PUBLIC   lim_hdf     ! called by lim_trp 
    2929 
    30    LOGICAL  ::   linit = .TRUE.              ! initialization flag (set to flase after the 1st call) 
    31    REAL(wp) ::   epsi04 = 1.e-04              ! constant 
     30   LOGICAL  ::   linit = .TRUE.                             ! initialization flag (set to flase after the 1st call) 
    3231   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   efact   ! metric coefficient 
    3332 
     
    5453      REAL(wp), DIMENSION(jpi,jpj), INTENT( inout ) ::   ptab    ! Field on which the diffusion is applied 
    5554      ! 
    56       INTEGER                           ::  ji, jj                   ! dummy loop indices 
    57       INTEGER                           ::  its, iter, ierr          ! local integers 
    58       REAL(wp)                          ::   zalfa, zrlxint, zconv   ! local scalars 
    59       REAL(wp), POINTER, DIMENSION(:,:) ::   zrlx, zflu, zflv, zdiv0, zdiv, ztab0 
    60       CHARACTER(lc)                     ::   charout                 ! local character 
    61       REAL(wp), PARAMETER               :: zrelax = 0.5_wp           ! relaxation constant for iterative procedure 
     55      INTEGER                           ::  ji, jj                    ! dummy loop indices 
     56      INTEGER                           ::  iter, ierr           ! local integers 
     57      REAL(wp)                          ::  zrlxint, zconv     ! local scalars 
     58      REAL(wp), POINTER, DIMENSION(:,:) ::  zrlx, zflu, zflv, zdiv0, zdiv, ztab0 
     59      CHARACTER(lc)                     ::  charout                   ! local character 
     60      REAL(wp), PARAMETER               ::  zrelax = 0.5_wp           ! relaxation constant for iterative procedure 
     61      REAL(wp), PARAMETER               ::  zalfa  = 0.5_wp           ! =1.0/0.5/0.0 = implicit/Cranck-Nicholson/explicit 
     62      INTEGER , PARAMETER               ::  its    = 100              ! Maximum number of iteration 
    6263      !!------------------------------------------------------------------- 
    6364       
     
    7879      ENDIF 
    7980      !                             ! Time integration parameters 
    80       zalfa = 0.5_wp                      ! =1.0/0.5/0.0 = implicit/Cranck-Nicholson/explicit 
    81       its   = 100                         ! Maximum number of iteration 
    8281      ! 
    8382      ztab0(:, : ) = ptab(:,:)      ! Arrays initialization 
     
    9291      iter  = 0 
    9392      ! 
    94       DO WHILE( zconv > ( 2._wp * epsi04 ) .AND. iter <= its )   ! Sub-time step loop 
     93      DO WHILE( zconv > ( 2._wp * 1.e-04 ) .AND. iter <= its )   ! Sub-time step loop 
    9594         ! 
    9695         iter = iter + 1                                 ! incrementation of the sub-time step number 
  • 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 
  • branches/2015/dev_r5044_CNRS_LIM3CLEAN/NEMOGCM/NEMO/OPA_SRC/BDY/bdyice_lim.F90

    r5079 r5080  
    230230               e_s(ji,jj,jk,jl) = rswitch * rhosn * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus ) 
    231231               ! Multiply by volume, so that heat content in J/m2 
    232                e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) / nlay_s 
     232               e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s 
    233233            END DO 
    234234            DO jk = 1, nlay_i 
     
    240240                  - rcp      * ( ztmelts - rt0 ) ) 
    241241               ! Mutliply by ice volume, and divide by number of layers to get heat content in J/m2 
    242                e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) / nlay_i 
     242               e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * a_i(ji,jj,jl) * ht_i(ji,jj,jl) * r1_nlay_i 
    243243            END DO 
    244244 
Note: See TracChangeset for help on using the changeset viewer.