Ignore:
Timestamp:
2013-11-20T09:48:15+01:00 (7 years ago)
Author:
clem
Message:

update lim3 for bdy purpose

Location:
branches/2013/dev_LOCEAN_CMCC_INGV_MERC_UKMO_2013/NEMOGCM/NEMO/LIM_SRC_3
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • branches/2013/dev_LOCEAN_CMCC_INGV_MERC_UKMO_2013/NEMOGCM/NEMO/LIM_SRC_3/limcat_1D.F90

    r4161 r4269  
    4545      REAL(wp), DIMENSION(:,:), INTENT(inout) ::   zht_i, zht_s, za_i ! output ice/snow variables 
    4646      REAL(wp) :: epsi06 = 1.0e-6 
    47    
     47      REAL(wp) ::   zc1, zc2, zc3, zx1   ! local scalars 
     48      REAL(wp), DIMENSION(0:jpl)   ::   zhi_max         !:Boundary of ice thickness categories in thickness space 
     49  
    4850      IF( nn_timing == 1 )  CALL timing_start('limcat_1D') 
    4951      !-------------------------------------------------------------------- 
     
    6668      !         fulfills ice volume concervation between input and output (ztests=4)  
    6769      !-------------------------------------------------------------------------------------- 
     70 
     71      !- Thickness categories boundaries 
     72      !    hi_max is calculated in iceini.F90 but since limcat_1D.F90 routine  
     73      !    is called before (in bdydta.F90), one must recalculate it.    
     74      !    Note clem: there may be a way of doing things cleaner 
     75      !---------------------------------- 
     76      zhi_max(:) = 0._wp 
     77      zc1 =  3._wp / REAL( jpl , wp ) ; zc2 = 10._wp * zc1 ; zc3 =  3._wp 
     78      DO jl = 1, jpl 
     79         zx1 = REAL( jl-1 , wp ) / REAL( jpl , wp ) 
     80         zhi_max(jl) = zhi_max(jl-1) + zc1 + zc2 * ( 1._wp + TANH( zc3 * ( zx1 - 1._wp ) ) ) 
     81      END DO 
    6882  
    6983      ! ---------------------------------------- 
     
    100114               ! Fill ice thicknesses except the last one (i_fill) by (hmax-hmin)/2  
    101115               DO jl = 1, i_fill - 1 
    102                   zht_i(ji,jl) = ( hi_max(jl) + hi_max(jl-1) ) / 2. 
     116                  zht_i(ji,jl) = ( zhi_max(jl) + zhi_max(jl-1) ) / 2. 
    103117               END DO 
    104118                
     
    106120               jl0 = i_fill 
    107121               DO jl = 1, i_fill 
    108                   IF ( ( zhti(ji) >= hi_max(jl-1) ) .AND. ( zhti(ji) < hi_max(jl) ) ) THEN 
     122                  IF ( ( zhti(ji) >= zhi_max(jl-1) ) .AND. ( zhti(ji) < zhi_max(jl) ) ) THEN 
    109123                     jl0 = jl 
    110124           CYCLE 
     
    141155             
    142156            ! Test 3: thickness of the last category is in-bounds ? 
    143             IF ( zht_i(ji,i_fill) >= hi_max(i_fill-1) ) ztest_3 = 1 
     157            IF ( zht_i(ji,i_fill) >= zhi_max(i_fill-1) ) ztest_3 = 1 
    144158             
    145159            ! Test 4: positivity of ice concentrations 
  • branches/2013/dev_LOCEAN_CMCC_INGV_MERC_UKMO_2013/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90

    r4205 r4269  
    426426               !-Calculate stress tensor components zs1 and zs2  
    427427               !-at centre of grid cells (see section 3.5 of CICE user's guide). 
    428                zs1(ji,jj) = ( zs1(ji,jj) - dtotel*( ( 1._wp - alphaevp) * zs1(ji,jj) +   & 
    429                   &          ( delta / deltat(ji,jj) - zdd(ji,jj) / deltat(ji,jj) ) * zpresh(ji,jj) ) )  &        
    430                   &          / ( 1._wp + alphaevp * dtotel ) 
    431  
    432                zs2(ji,jj) = ( zs2(ji,jj) - dtotel * ( ( 1._wp - alphaevp ) * ecc2 * zs2(ji,jj) -   & 
    433                              zdt(ji,jj) / deltat(ji,jj) * zpresh(ji,jj) ) )   & 
    434                   &          / ( 1._wp + alphaevp * ecc2 * dtotel ) 
     428               !zs1(ji,jj) = ( zs1(ji,jj) - dtotel*( ( 1._wp - alphaevp) * zs1(ji,jj) +   & 
     429               !   &          ( delta / deltat(ji,jj) - zdd(ji,jj) / deltat(ji,jj) ) * zpresh(ji,jj) ) )  &        
     430               !   &          / ( 1._wp + alphaevp * dtotel ) 
     431 
     432               !zs2(ji,jj) = ( zs2(ji,jj) - dtotel * ( ( 1._wp - alphaevp ) * ecc2 * zs2(ji,jj) -   & 
     433               !              zdt(ji,jj) / deltat(ji,jj) * zpresh(ji,jj) ) )   & 
     434               !   &          / ( 1._wp + alphaevp * ecc2 * dtotel ) 
    435435 
    436436               ! new formulation from S. Bouillon to help stabilizing the code (no need of alphaevp) 
    437                !zs1(ji,jj) = ( zs1(ji,jj) + dtotel * ( ( zdd(ji,jj) / deltat(ji,jj) - delta / deltat(ji,jj) )  & 
    438                !   &         * zpresh(ji,jj) ) ) / ( 1._wp + dtotel ) 
    439                !zs2(ji,jj) = ( zs2(ji,jj) + dtotel * ( ecci * zdt(ji,jj) / deltat(ji,jj) * zpresh(ji,jj) ) )  & 
    440                !   &         / ( 1._wp + dtotel ) 
     437               zs1(ji,jj) = ( zs1(ji,jj) + dtotel * ( ( zdd(ji,jj) / deltat(ji,jj) - delta / deltat(ji,jj) )  & 
     438                  &         * zpresh(ji,jj) ) ) / ( 1._wp + dtotel ) 
     439               zs2(ji,jj) = ( zs2(ji,jj) + dtotel * ( ecci * zdt(ji,jj) / deltat(ji,jj) * zpresh(ji,jj) ) )  & 
     440                  &         / ( 1._wp + dtotel ) 
    441441 
    442442            END DO 
     
    472472 
    473473               !-Calculate stress tensor component zs12 at corners (see section 3.5 of CICE user's guide). 
    474                zs12(ji,jj) = ( zs12(ji,jj) - dtotel * ( (1.0-alphaevp) * ecc2 * zs12(ji,jj) - zds(ji,jj) /  & 
    475                   &          ( 2._wp * deltac(ji,jj) ) * zpreshc(ji,jj) ) )  & 
    476                   &          / ( 1._wp + alphaevp * ecc2 * dtotel )  
     474               !zs12(ji,jj) = ( zs12(ji,jj) - dtotel * ( (1.0-alphaevp) * ecc2 * zs12(ji,jj) - zds(ji,jj) /  & 
     475               !   &          ( 2._wp * deltac(ji,jj) ) * zpreshc(ji,jj) ) )  & 
     476               !   &          / ( 1._wp + alphaevp * ecc2 * dtotel )  
    477477 
    478478               ! new formulation from S. Bouillon to help stabilizing the code (no need of alphaevp) 
    479                !zs12(ji,jj) = ( zs12(ji,jj) + dtotel *  & 
    480                !   &          ( ecci * zds(ji,jj) / ( 2._wp * deltac(ji,jj) ) * zpreshc(ji,jj) ) )  & 
    481                !   &          / ( 1.0 + dtotel )  
     479               zs12(ji,jj) = ( zs12(ji,jj) + dtotel *  & 
     480                  &          ( ecci * zds(ji,jj) / ( 2._wp * deltac(ji,jj) ) * zpreshc(ji,jj) ) )  & 
     481                  &          / ( 1.0 + dtotel )  
    482482 
    483483            END DO ! ji 
Note: See TracChangeset for help on using the changeset viewer.