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 827 for branches – NEMO

Changeset 827 for branches


Ignore:
Timestamp:
2008-03-06T15:39:49+01:00 (16 years ago)
Author:
ctlod
Message:

dev_002_LIM : adapt ocean modules to both LIM 2.0 & 3.0 using key_lim2 & key_lim3

Location:
branches/dev_002_LIM/NEMO
Files:
2 deleted
8 edited

Legend:

Unmodified
Added
Removed
  • branches/dev_002_LIM/NEMO/LIM_SRC_2/limflx_2.F90

    r823 r827  
    1717   USE flx_oce          ! sea-ice/ocean forcings variables 
    1818   USE ice_2             ! LIM sea-ice variables 
    19    USE flxblk_2           ! bulk formulea 
     19   USE flxblk             ! bulk formulea 
    2020   USE lbclnk           ! ocean lateral boundary condition 
    2121   USE in_out_manager   ! I/O manager 
    22    USE albedo_2         ! albedo parameters 
     22   USE albedo           ! albedo parameters 
    2323   USE prtctl           ! Print control 
    2424 
     
    193193      !  2) Computation of snow/ice and ocean albedo   ! 
    194194      !------------------------------------------------! 
    195       CALL flx_blk_albedo_2( zalb, zalcn, zalbp, zaldum ) 
     195      CALL flx_blk_albedo( zalb, zalcn, zalbp, zaldum ) 
    196196 
    197197      alb_ice(:,:) =  0.5 * zalbp(:,:) + 0.5 * zalb (:,:)   ! Ice albedo                        
  • branches/dev_002_LIM/NEMO/OPA_SRC/SBC/albedo.F90

    r826 r827  
    3737      c1     = 0.05  ,     &   ! constants values 
    3838      c2     = 0.10  ,     & 
     39#if defined key_lim3 
     40      albice = 0.53  ,     &   !  albedo of melting ice in the arctic and antarctic (Shine & Hendersson-Sellers) 
     41#elif defined key_lim2 
    3942      albice = 0.50  ,     &   !  albedo of melting ice in the arctic and antarctic (Shine & Hendersson-Sellers) 
     43#endif 
    4044      cgren  = 0.06  ,     &   !  correction of the snow or ice albedo to take into account 
    4145                               !  effects of cloudiness (Grenfell & Perovich, 1984) 
     
    5357CONTAINS 
    5458 
    55 #if defined key_lim3 
    56    !!---------------------------------------------------------------------- 
    57    !!   'key_lim3'                                         LIM ice model 
     59#if defined key_lim3 || defined key_lim2 
     60   !!---------------------------------------------------------------------- 
     61   !!   'key_lim3' OR 'key_lim2'               LIM 2.0 or LIM 3.0 ice model 
    5862   !!---------------------------------------------------------------------- 
    5963 
     
    7579      !!  8.0   !  01-04  (LIM 1.0) 
    7680      !!  8.5   !  03-07  (C. Ethe, G. Madec)  Optimization (old name:shine) 
     81      !!  9.0   !  01-06  (M. Vancoppenolle) LIM 3.0 
    7782      !!---------------------------------------------------------------------- 
    7883      !! * Modules used 
     84#if defined key_lim3 
    7985      USE ice                   ! ??? 
     86#elif defined key_lim2 
     87      USE ice_2                 ! ??? 
     88#endif 
    8089 
    8190      !! * Arguments 
     91#if defined key_lim3 
     92      REAL(wp), DIMENSION(jpi,jpj,jpl), INTENT(out) ::  & 
     93#elif defined key_lim2 
    8294      REAL(wp), DIMENSION(jpi,jpj), INTENT(out) ::  & 
     95#endif 
    8396         palb         ,     &    !  albedo of ice under overcast sky 
     97         palbp                   !  albedo of ice under clear sky  
     98      REAL(wp), DIMENSION(jpi,jpj), INTENT(out) ::  & 
    8499         palcn        ,     &    !  albedo of ocean under overcast sky 
    85          palbp        ,     &    !  albedo of ice under clear sky  
    86100         palcnp                  !  albedo of ocean under clear sky 
    87101 
    88102      !! * Local variables 
    89103      INTEGER ::    & 
    90          ji, jj                   ! dummy loop indices 
     104         ji, jj, jl               ! dummy loop indices 
    91105      REAL(wp) ::   &  
    92106         zmue14         ,     &   !  zmue**1.4 
     
    96110         zalbpic        ,     &   !  albedo of snow/ice system when ice is free of snow 
    97111         zithsn         ,     &   !  = 1 for hsn >= 0 ( ice is cov. by snow ) ; = 0 otherwise (ice is free of snow) 
    98          zitmlsn        ,     &   !  = 1 freezinz snow (sist >=rt0_snow) ; = 0 melting snow (sist<rt0_snow) 
     112         zitmlsn        ,     &   !  = 1 freezinz snow (t_su >=rt0_snow) ; = 0 melting snow (t_su<rt0_snow) 
    99113         zihsc1         ,     &   !  = 1 hsn <= c1 ; = 0 hsn > c1 
    100114         zihsc2                   !  = 1 hsn >= c2 ; = 0 hsn < c2 
     115#if defined key_lim3 
     116      REAL(wp), DIMENSION(jpi,jpj,jpl) ::  & 
     117#elif defined key_lim2 
    101118      REAL(wp), DIMENSION(jpi,jpj) ::  & 
     119#endif 
    102120         zalbfz         ,     &   !  ( = alphdi for freezing ice ; = albice for melting ice ) 
    103121         zficeth                  !  function of ice thickness 
     122#if defined key_lim3 
     123      LOGICAL , DIMENSION(jpi,jpj,jpl) ::  & 
     124#elif defined key_lim2 
    104125      LOGICAL , DIMENSION(jpi,jpj) ::  & 
     126#endif 
    105127         llmask 
    106128      !!--------------------------------------------------------------------- 
     
    112134      !  Computation of  zficeth 
    113135      !--------------------------  
    114        
     136#if defined key_lim3 
     137      llmask = (ht_s(:,:,:) == 0.e0) .AND. ( t_su(:,:,:) >= rt0_ice ) 
     138#elif defined key_lim2       
    115139      llmask = (hsnif == 0.e0) .AND. ( sist >= rt0_ice ) 
     140#endif 
    116141      WHERE ( llmask )   !  ice free of snow and melts 
    117142         zalbfz = albice 
     
    120145      END WHERE 
    121146       
     147#if defined key_lim3 
     148      DO jl = 1, jpl 
     149         DO jj = 1, jpj 
     150            DO ji = 1, jpi 
     151               IF( ht_i(ji,jj,jl) > 1.5 ) THEN 
     152                  zficeth(ji,jj,jl) = zalbfz(ji,jj,jl) 
     153               ELSEIF( ht_i(ji,jj,jl) > 1.0  .AND. ht_i(ji,jj,jl) <= 1.5 ) THEN 
     154                  zficeth(ji,jj,jl) = 0.472 + 2.0 * ( zalbfz(ji,jj,jl) - 0.472 ) * ( ht_i(ji,jj,jl) - 1.0 ) 
     155               ELSEIF( ht_i(ji,jj,jl) > 0.05 .AND. ht_i(ji,jj,jl) <= 1.0 ) THEN 
     156                  zficeth(ji,jj,jl) = 0.2467 + 0.7049 * ht_i(ji,jj,jl)                               & 
     157                     &                    - 0.8608 * ht_i(ji,jj,jl) * ht_i(ji,jj,jl)                 & 
     158                     &                    + 0.3812 * ht_i(ji,jj,jl) * ht_i(ji,jj,jl) * ht_i (ji,jj,jl) 
     159               ELSE 
     160                  zficeth(ji,jj,jl) = 0.1 + 3.6 * ht_i(ji,jj,jl)  
     161               ENDIF 
     162            END DO 
     163         END DO 
     164      END DO 
     165#elif defined key_lim2       
    122166      DO jj = 1, jpj 
    123167         DO ji = 1, jpi 
     
    135179         END DO 
    136180      END DO 
     181#endif 
    137182       
    138183      !-----------------------------------------------  
     
    142187      !    Albedo of snow-ice for clear sky. 
    143188      !-----------------------------------------------     
     189#if defined key_lim3 
     190      DO jl = 1, jpl 
     191         DO jj = 1, jpj 
     192            DO ji = 1, jpi 
     193               !  Case of ice covered by snow.              
     194             
     195               !  freezing snow         
     196               zihsc1       = 1.0 - MAX ( zzero , SIGN ( zone , - ( ht_s(ji,jj,jl) - c1 ) ) ) 
     197               zalbpsnf     = ( 1.0 - zihsc1 ) * ( zficeth(ji,jj,jl) + ht_s(ji,jj,jl) * ( alphd - zficeth(ji,jj,jl) ) / c1 ) & 
     198                  &                 + zihsc1   * alphd   
     199 
     200               !  melting snow                 
     201               zihsc2       = MAX ( zzero , SIGN ( zone , ht_s(ji,jj,jl) - c2 ) ) 
     202               zalbpsnm     = ( 1.0 - zihsc2 ) * ( albice + ht_s(ji,jj,jl) * ( alphc - albice ) / c2 )                 & 
     203                  &                 + zihsc2   * alphc  
     204 
     205 
     206               zitmlsn      =  MAX ( zzero , SIGN ( zone , t_su(ji,jj,jl) - rt0_snow ) )    
     207               zalbpsn      =  zitmlsn * zalbpsnm + ( 1.0 - zitmlsn ) * zalbpsnf 
     208             
     209               !  Case of ice free of snow. 
     210               zalbpic      = zficeth(ji,jj,jl)  
     211             
     212               ! albedo of the system    
     213               zithsn       = 1.0 - MAX ( zzero , SIGN ( zone , - ht_s(ji,jj,jl) ) ) 
     214               palbp(ji,jj,jl) =  zithsn * zalbpsn + ( 1.0 - zithsn ) *  zalbpic 
     215            END DO 
     216         END DO 
     217      END DO 
     218       
     219      !    Albedo of snow-ice for overcast sky. 
     220      !----------------------------------------------   
     221      palb(:,:,:)   = palbp(:,:,:) + cgren       ! Oberhuber correction 
     222 
     223#elif defined key_lim2       
     224 
    144225      DO jj = 1, jpj 
    145226         DO ji = 1, jpi 
     
    170251      !----------------------------------------------   
    171252      palb(:,:)   = palbp(:,:) + cgren                                            
     253#endif 
    172254       
    173255      !-------------------------------------------- 
  • branches/dev_002_LIM/NEMO/OPA_SRC/SBC/bulk.F90

    r826 r827  
    1616   USE ocfzpt          ! ocean freezing point 
    1717   USE flxblk          ! bulk formulae 
    18    USE flxblk_2        ! bulk formulae 
    1918   USE blk_oce         ! bulk variable  
    2019   USE flx_oce 
     
    8887 
    8988         zsst(:,:) = gsst(:,:) / REAL( nfbulk ) * tmask(:,:,1) 
    90          IF ( .NOT. lk_lim2 )  CALL flx_blk( zsst )     
    9189 
    92 #if defined key_lim2 
    93          CALL flx_blk_2( zsst )     
    94 #endif 
     90         CALL flx_blk( zsst )     
    9591      
    9692         gsst(:,:) = 0.     
  • branches/dev_002_LIM/NEMO/OPA_SRC/SBC/flx_oce.F90

    r826 r827  
    1313   !! * Modules used 
    1414   USE par_oce          ! ocean parameters 
     15# if defined key_lim3 
     16   USE par_ice 
     17# endif 
    1518 
    1619   IMPLICIT NONE 
     
    3134 
    3235#elif defined key_lim3 || defined key_lim2 || defined key_flx_bulk_monthly || defined key_flx_bulk_daily || defined key_flx_core 
     36 
     37#if defined key_lim3  
     38   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpl)    ::   &  !: 
     39#elif defined key_lim2  
    3340   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)    ::   &  !: 
     41#endif 
    3442      qsr_ice  ,      &  !: solar flux over ice 
     43      tn_ice   ,      &  !: ice surface temperature 
     44      qnsr_ice ,      &  !: total non solar heat flux (Longwave downward radiation) over ice 
     45      dqns_ice           !: total non solar sensibility over ice (LW+SEN+LA) 
     46   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)    ::   &  !: 
    3547      qsr_oce  ,      &  !: solar flux over ocean 
    3648      qnsr_oce ,      &  !: total non solar heat flux (Longwave downward radiation) over ocean  
    37       qnsr_ice ,      &  !: total non solar heat flux (Longwave downward radiation) over ice 
    3849      tprecip  ,      &  !: total precipitation ( or liquid precip minus evaporation in coupled mode) 
    3950      sprecip  ,      &  !: solid (snow) precipitation 
    40       dqns_ice ,      &  !: total non solar sensibility over ice (LW+SEN+LA) 
    41       tn_ice   ,      &  !: ice surface temperature 
    4251      evap     ,      &  !: evaporation over ocean 
    4352      fr1_i0   ,      &  !: 1st part of the fraction of sol. rad.  which penetrate inside the ice cover 
    44       fr2_i0   ,      &  !: 2nd part of the fraction of sol. rad.  which penetrate inside the ice cover  
     53      fr2_i0             !: 2nd part of the fraction of sol. rad.  which penetrate inside the ice cover  
    4554#if ! defined key_coupled 
     55#if defined key_lim3  
     56   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpl)    ::   &  !: 
     57#elif defined key_lim2  
     58   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)    ::   &  !: 
     59#endif 
    4660      qla_ice  ,      &  !: latent flux over ice   
    4761      dqla_ice           !: latent sensibility over ice 
    48 #else 
     62#elif  
     63   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)    ::   &  !: 
    4964      rrunoff  ,      &  !: runoff 
    5065      calving  ,      &  !: calving 
  • branches/dev_002_LIM/NEMO/OPA_SRC/SBC/flxblk.F90

    r826 r827  
    2525   USE albedo 
    2626   USE prtctl          ! Print control 
    27  
     27#if defined key_lim3 
     28   USE ice 
     29#elif defined key_lim2 
     30   USE ice_2 
     31#endif 
    2832   IMPLICIT NONE 
    2933   PRIVATE 
     
    7983 
    8084CONTAINS 
     85#if ! defined key_lim2 
     86 
     87   SUBROUTINE flx_blk( psst ) 
     88      !!--------------------------------------------------------------------------- 
     89      !!                     ***  ROUTINE flx_blk  *** 
     90      !!                  
     91      !!  ** Purpose :   Computation of the heat fluxes at ocean and snow/ice 
     92      !!       surface the solar heat at ocean and snow/ice surfaces and the  
     93      !!       sensitivity of total heat fluxes to the SST variations 
     94      !!          
     95      !!  ** Method  :   The flux of heat at the ice and ocean surfaces are derived 
     96      !!       from semi-empirical ( or bulk ) formulae which relate the flux to  
     97      !!       the properties of the surface and of the lower atmosphere. Here, we 
     98      !!       follow the work of Oberhuber, 1988    
     99      !! 
     100      !!  ** Action  :   call flx_blk_albedo to compute ocean and ice albedo  
     101      !!          computation of snow precipitation 
     102      !!          computation of solar flux at the ocean and ice surfaces 
     103      !!          computation of the long-wave radiation for the ocean and sea/ice 
     104      !!          computation of turbulent heat fluxes over water and ice 
     105      !!          computation of evaporation over water 
     106      !!          computation of total heat fluxes sensitivity over ice (dQ/dT) 
     107      !!          computation of latent heat flux sensitivity over ice (dQla/dT) 
     108      !! 
     109      !! History : 
     110      !!   8.0  !  97-06  (Louvain-La-Neuve)  Original code 
     111      !!   8.5  !  02-09  (C. Ethe , G. Madec )  F90: Free form and module 
     112      !!---------------------------------------------------------------------- 
     113      !! * Arguments 
     114      REAL(wp), INTENT( in ), DIMENSION(jpi,jpj) ::   & 
     115         &                          psst      ! Sea Surface Temperature  
     116 
     117      !! * Local variables 
     118      INTEGER  ::             & 
     119         ji, jj, jl, jt    ,  &  ! dummy loop indices 
     120         indaet            ,  &  !  = -1, 0, 1 for odd, normal and leap years resp. 
     121         iday              ,  &  ! integer part of day 
     122         indxb             ,  &  ! index for budyko coefficient 
     123         indxc                   ! index for cloud depth coefficient 
     124 
     125      REAL(wp)  ::            &  
     126         zalat , zclat     ,  &  ! latitude in degrees  
     127         zmt1, zmt2, zmt3  ,  &  ! tempory air temperatures variables 
     128         ztatm3, ztatm4    ,  &  ! power 3 and 4 of air temperature 
     129         z4tatm3           ,  &  ! 4 * ztatm3 
     130         zcmue             ,  &  ! cosine of local solar altitude 
     131         zcmue2            ,  &  ! root of zcmue1  
     132         zscmue            ,  &  ! square-root of zcmue1  
     133         zpcmue            ,  &  ! zcmue1**1.4 
     134         zdecl             ,  &  ! solar declination 
     135         zsdecl , zcdecl   ,  &  ! sine and cosine of solar declination  
     136         zalbo             ,  &  ! albedo of sea-water 
     137         zalbi             ,  &  ! albedo of ice 
     138         ztamr             ,  &  ! air temperature minus triple point of water (rtt) 
     139         ztaevbk           ,  &  ! part of net longwave radiation 
     140         zevi , zevo       ,  &  ! vapour pressure of ice and ocean  
     141         zind1,zind2,zind3 ,  &  ! switch for testing the values of air temperature 
     142         zinda             ,  &  ! switch for testing the values of sea ice cover 
     143         zpis2             ,  &  ! pi / 2 
     144         z2pi                    ! 2 * pi  
     145 
     146      REAL(wp)  ::            &  
     147         zxday             ,  &  ! day of year 
     148         zdist             ,  &  ! distance between the sun and the earth during the year 
     149         zdaycor           ,  &  ! corr. factor to take into account the variation of  
     150         !                       ! zday when calc. the solar rad.     
     151         zesi, zeso        ,  &  ! vapour pressure of ice and ocean at saturation 
     152         zesi2             ,  &  ! root of zesi  
     153         zqsato            ,  &  ! humidity close to the ocean surface (at saturation)    
     154         zqsati            ,  &  ! humidity close to the ice surface (at saturation)  
     155         zqsati2           ,  &  ! root of  zqsati  
     156         zdesidt           ,  &  ! derivative of zesi, function of ice temperature 
     157         zdteta            ,  &  ! diff. betw. sst and air temperature 
     158         zdeltaq           ,  &  ! diff. betw. spec. hum. and hum. close to the surface 
     159         ztvmoy, zobouks   ,  &  ! tempory scalars 
     160         zpsims, zpsihs, zpsils, zobouku, zxins, zpsimu ,  &  
     161         zpsihu, zpsilu, zstab,zpsim, zpsih, zpsil      ,  &  
     162         zvatmg, zcmn, zchn, zcln, zcmcmn, zdenum       ,  &  
     163         zdtetar, ztvmoyr, zlxins, zcmn2, zchcm, zclcm , zcoef 
     164 
     165      REAL(wp)  ::            &  
     166         zrhova            ,  &  ! air density per wind speed 
     167         zcsho , zcleo     ,  &  ! transfer coefficient over ocean 
     168         zcshi , zclei     ,  &  ! transfer coefficient over ice-free 
     169         zrhovacleo        ,  &  ! air density per wind speed per transfer coef. 
     170         zrhovacsho, zrhovaclei, zrhovacshi, &  
     171         ztice3            ,  &  ! power 3 of ice temperature 
     172         zticemb, zticemb2 ,  &  ! tempory air temperatures variables 
     173         zdqlw_ice         ,  &  ! sensitivity of long-wave flux over ice 
     174         zdqsb_ice         ,  &  ! sensitivity of sensible heat flux over ice 
     175         zdqla_ice         ,  &  ! sensitivity of latent heat flux over ice 
     176         zdl, zdr                ! fractionnal part of latitude 
     177      REAL(wp), DIMENSION(jpi,jpj) :: &  
     178         zpatm            ,  &   ! atmospheric pressure 
     179         zqatm            ,  &   ! specific humidity 
     180         zes              ,  &   ! vapour pressure at saturation 
     181         zev, zevsqr      ,  &   ! vapour pressure and his square-root 
     182         zrhoa            ,  &   ! air density 
     183         ztatm            ,  &   ! air temperature in Kelvins 
     184         zfrld            ,  &   ! fraction of sea ice cover  
     185         zcatm1           ,  &   ! fraction of cloud 
     186         zcldeff                 ! correction factor to account cloud effect 
     187      REAL(wp), DIMENSION(jpi,jpj) ::   &  
     188         zalbocsd         ,  &   ! albedo of ocean 
     189         zalboos          ,  &   ! albedo of ocean under overcast sky 
     190         zalbomu          ,  &   ! albedo of ocean when zcmue is 0.4 
     191         zqsro            ,  &   ! solar radiation over ocean 
     192         zqsrics          ,  &   ! solar radiation over ice under clear sky 
     193         zqsrios          ,  &   ! solar radiation over ice under overcast sky 
     194         zcldcor          ,  &   ! cloud correction 
     195         zlsrise, zlsset  ,  &   ! sunrise and sunset 
     196         zlmunoon         ,  &   ! local noon solar altitude 
     197         zdlha            ,  &   ! length of the ninstr segments of the solar day 
     198         zps              ,  &   ! sine of latitude per sine of solar decli. 
     199         zpc                     ! cosine of latitude per cosine of solar decli.  
     200      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   &  
     201         zalbics          ,  &   ! albedo of ice under clear sky 
     202         zalbios                 ! albedo of ice under overcast sky 
     203 
     204      REAL(wp), DIMENSION(jpi,jpj) ::   &  
     205         zqlw_oce         ,  &   ! long-wave heat flux over ocean 
     206         zqla_oce         ,  &   ! latent heat flux over ocean 
     207         zqsb_oce                ! sensible heat flux over ocean 
     208  
     209      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   &  
     210         zqlw_ice         ,  &   ! long-wave heat flux over ice 
     211         zqla_ice         ,  &   ! latent heat flux over ice 
     212         zqsb_ice                ! sensible heat flux over ice 
     213  
     214      REAL(wp), DIMENSION(jpi,jpj,jpintsr) ::    & 
     215         zlha             ,  &   ! local hour angle 
     216         zalbocs          ,  &   ! tempory var. of ocean albedo under clear sky 
     217         zsqsro           ,  &   ! tempory var. of solar rad. over ocean  
     218         zsqsrics         ,  &   ! temp. var. of solar rad. over ice under clear sky 
     219         zsqsrios                ! temp. var. of solar rad. over ice under overcast sky 
     220      !!--------------------------------------------------------------------- 
     221 
     222      !--------------------- 
     223      !   Initilization    ! 
     224      !--------------------- 
     225#if ! defined key_lim3 
     226      tn_ice(:,:) = psst(:,:) 
     227#endif 
     228 
     229      !  Determine cloud optical depths as a function of latitude (Chou et al., 1981). 
     230      !  and the correction factor for taking into account  the effect of clouds  
     231      !------------------------------------------------------ 
     232      IF( lbulk_init ) THEN 
     233         DO jj = 1, jpj   
     234            DO ji = 1 , jpi 
     235               zalat          = ( 90.e0 - ABS( gphit(ji,jj) ) ) /  5.e0 
     236               zclat          = ( 95.e0 -      gphit(ji,jj)   ) / 10.e0 
     237               indxb          = 1 + INT( zalat )  
     238               !  correction factor to account for the effect of clouds  
     239               sbudyko(ji,jj) = budyko(indxb)   
     240               indxc          = 1 + INT( zclat )   
     241               zdl            = zclat - INT( zclat )  
     242               zdr            = 1.0 - zdl 
     243               stauc(ji,jj)   = zdr * tauco( indxc ) + zdl * tauco( indxc + 1 )  
     244            END DO 
     245         END DO 
     246         IF( nleapy == 1 ) THEN 
     247            yearday = 366.e0 
     248         ELSE IF( nleapy == 0 ) THEN 
     249            yearday = 365.e0 
     250         ELSEIF( nleapy == 30) THEN 
     251            yearday = 360.e0 
     252         ENDIF 
     253         lbulk_init = .FALSE. 
     254      ENDIF 
     255 
     256      zqlw_oce(:,:) = 0.e0 
     257      zqla_oce(:,:) = 0.e0 
     258      zqsb_oce(:,:) = 0.e0 
     259      zqlw_ice(:,:,:) = 0.e0 
     260      zqla_ice(:,:,:) = 0.e0 
     261      zqsb_ice(:,:,:) = 0.e0 
     262 
     263      zpis2       = rpi / 2. 
     264      z2pi        = 2. * rpi 
     265 
     266 !CDIR NOVERRCHK 
     267      DO jj = 1, jpj 
     268 !CDIR NOVERRCHK 
     269         DO ji = 1, jpi 
     270 
     271            ztatm (ji,jj) = 273.15 + tatm  (ji,jj)  !  air temperature in Kelvins  
     272            zcatm1(ji,jj) = 1.0    - catm  (ji,jj)  !  fractional cloud cover 
     273            zfrld (ji,jj) = 1.0    - freeze(ji,jj)  !  fractional sea ice cover 
     274            zpatm(ji,jj)  = 101000.               !  pressure  
     275       
     276            !  Computation of air density, obtained from the equation of state for dry air.  
     277            zrhoa(ji,jj) = zpatm(ji,jj) / ( 287.04 * ztatm(ji,jj) ) 
     278       
     279            !  zes : Saturation water vapour 
     280            ztamr = ztatm(ji,jj) - rtt 
     281            zmt1  = SIGN( 17.269, ztamr ) 
     282            zmt2  = SIGN( 21.875, ztamr ) 
     283            zmt3  = SIGN( 28.200, -ztamr ) 
     284            zes(ji,jj) = 611.0 * EXP (  ABS( ztamr ) * MIN ( zmt1, zmt2 )   & 
     285               &                      / ( ztatm(ji,jj) - 35.86  + MAX( zzero, zmt3 ) ) ) 
     286 
     287            !  zev : vapour pressure  (hatm is relative humidity)   
     288            zev(ji,jj)   = hatm(ji,jj) * zes(ji,jj)  
     289            !  square-root of vapour pressure 
     290!CDIR NOVERRCHK 
     291            zevsqr(ji,jj) = SQRT( zev(ji,jj) * 0.01 ) 
     292            !  zqapb  : specific humidity  
     293            zqatm(ji,jj) = 0.622 * zev(ji,jj) / ( zpatm(ji,jj) - 0.378 * zev(ji,jj) ) 
     294 
     295 
     296            !---------------------------------------------------- 
     297            !   Computation of snow precipitation (Ledley, 1985) | 
     298            !---------------------------------------------------- 
     299 
     300            zmt1  =   253.0 - ztatm(ji,jj) 
     301            zmt2  = ( 272.0 - ztatm(ji,jj) ) / 38.0  
     302            zmt3  = ( 281.0 - ztatm(ji,jj) ) / 18.0 
     303            zind1 = MAX( zzero, SIGN( zone, zmt1 ) ) 
     304            zind2 = MAX( zzero, SIGN( zone, zmt2 ) ) 
     305            zind3 = MAX( zzero, SIGN( zone, zmt3 ) ) 
     306            ! total precipitation 
     307            tprecip(ji,jj) = watm(ji,jj) 
     308            ! solid  (snow) precipitation 
     309            sprecip(ji,jj) = tprecip(ji,jj) *       & 
     310               &             (           zind1      & 
     311               &               + ( 1.0 - zind1 ) * ( zind2 * ( 0.5 + zmt2 ) + ( 1.0 - zind2 ) *  zind3 * zmt3 ) )  
     312         END DO 
     313      END DO 
     314 
     315      !---------------------------------------------------------- 
     316      !   Computation of albedo (need to calculates heat fluxes)| 
     317      !----------------------------------------------------------- 
     318       
     319      CALL flx_blk_albedo( zalbios, zalboos, zalbics, zalbomu ) 
     320 
     321      !------------------------------------- 
     322      !   Computation of solar irradiance. | 
     323      !---------------------------------------- 
     324      indaet   = 1   
     325      !  compution of the day of the year at which the fluxes have to be calculate  
     326      !--The date corresponds to the middle of the time step. 
     327      zxday=nday_year + rdtbs2/rday 
     328 
     329      iday   = INT( zxday ) 
     330 
     331      IF(ln_ctl) CALL prt_ctl_info('declin : iday ', ivar1=iday, clinfo2=' nfbulk= ', ivar2=nfbulk) 
     332 
     333      !   computation of the solar declination, his sine and his cosine 
     334      CALL flx_blk_declin( indaet, iday, zdecl ) 
     335       
     336      zdecl    = zdecl * rad 
     337      zsdecl   = SIN( zdecl ) 
     338      zcdecl   = COS( zdecl ) 
     339       
     340      !  correction factor added for computation of shortwave flux to take into account the variation of 
     341      !  the distance between the sun and the earth during the year (Oberhuber 1988) 
     342      zdist    = zxday * z2pi / yearday 
     343      zdaycor  = 1.0 + 0.0013 * SIN( zdist ) + 0.0342 * COS( zdist ) 
     344 
     345!CDIR NOVERRCHK 
     346      DO jj = 1, jpj 
     347!CDIR NOVERRCHK 
     348         DO ji = 1, jpi 
     349            !  product of sine of latitude and sine of solar declination 
     350            zps     (ji,jj) = SIN( gphit(ji,jj) * rad ) * zsdecl 
     351            !  product of cosine of latitude and cosine of solar declination 
     352            zpc     (ji,jj) = COS( gphit(ji,jj) * rad ) * zcdecl 
     353            !  computation of the both local time of sunrise and sunset 
     354            zlsrise (ji,jj) = ACOS( - SIGN( zone, zps(ji,jj) ) * MIN( zone, SIGN( zone, zps(ji,jj) )  & 
     355               &                     * ( zps(ji,jj) / zpc(ji,jj) ) ) )  
     356            zlsset  (ji,jj) = - zlsrise(ji,jj) 
     357            !  dividing the solar day into jpintsr segments of length zdlha 
     358            zdlha   (ji,jj) = ( zlsrise(ji,jj) - zlsset(ji,jj) ) / REAL( jpintsr ) 
     359            !  computation of the local noon solar altitude 
     360            zlmunoon(ji,jj) = ASIN ( ( zps(ji,jj) + zpc(ji,jj) ) ) / rad 
     361             
     362            !  cloud correction taken from Reed (1977) (imposed lower than 1) 
     363            zcldcor (ji,jj) = MIN( zone, ( 1.e0 - 0.62 * catm(ji,jj) + 0.0019 * zlmunoon(ji,jj) ) ) 
     364         END DO 
     365      END DO 
     366 
     367         !  Computation of solar heat flux at each time of the day between sunrise and sunset.  
     368         !  We do this to a better optimisation of the code  
     369         !------------------------------------------------------        
     370      DO jl = 1, jpl 
     371 
     372!CDIR NOVERRCHK    
     373      DO jt = 1, jpintsr    
     374         zcoef = FLOAT( jt ) - 0.5 
     375!CDIR NOVERRCHK      
     376         DO jj = 1, jpj 
     377!CDIR NOVERRCHK 
     378            DO ji = 1, jpi 
     379               !  local hour angle 
     380               zlha (ji,jj,jt) = COS ( zlsrise(ji,jj) - zcoef * zdlha(ji,jj) ) 
     381 
     382               ! cosine of local solar altitude 
     383               zcmue              = MAX ( zzero ,   zps(ji,jj) + zpc(ji,jj) * zlha (ji,jj,jt)  ) 
     384               zcmue2             = 1368.0 * zcmue * zcmue 
     385               zscmue             = SQRT ( zcmue ) 
     386               zpcmue             = zcmue**1.4 
     387               ! computation of sea-water albedo (Payne, 1972) 
     388               zalbocs(ji,jj,jt)  = 0.05 / ( 1.1 * zpcmue + 0.15 ) 
     389               zalbo              = zcatm1(ji,jj) * zalbocs(ji,jj,jt) + catm(ji,jj) * zalboos(ji,jj) 
     390               ! solar heat flux absorbed at ocean surfaces (Zillman, 1972) 
     391               zevo               = zev(ji,jj) * 1.0e-05 
     392               zsqsro(ji,jj,jt)   =  ( 1.0 - zalbo ) * zdlha(ji,jj) * zcmue2                & 
     393                                   / ( ( zcmue + 2.7 ) * zevo + 1.085 * zcmue +  0.10 ) 
     394               !  solar heat flux absorbed at sea/ice surfaces  
     395               !  Formulation of Shine and Crane, 1984 adapted for high albedo surfaces  
     396 
     397               !  For clear sky         
     398               zevi               = zevo 
     399               zalbi              = zalbics(ji,jj,jl) 
     400               zsqsrics(ji,jj,jt) =  ( 1.0 - zalbi ) * zdlha(ji,jj) * zcmue2                & 
     401                  &                / ( ( 1.0 + zcmue ) * zevi + 1.2 * zcmue + 0.0455 ) 
     402 
     403               ! For overcast sky 
     404               zalbi              = zalbios(ji,jj,jl) 
     405               zsqsrios(ji,jj,jt) = zdlha(ji,jj) *                                                           & 
     406                  &                 ( ( 53.5 + 1274.5 * zcmue      ) *  zscmue * ( 1.0 - 0.996  * zalbi ) )  & 
     407                  &                 / (  1.0 + 0.139  * stauc(ji,jj) *           ( 1.0 - 0.9435 * zalbi ) ) 
     408            END DO 
     409         END DO 
     410      END DO 
     411 
     412 
     413      !  Computation of daily (between sunrise and sunset) solar heat flux absorbed  
     414      !  at the ocean and snow/ice surfaces. 
     415      !-------------------------------------------------------------------- 
     416 
     417      zalbocsd(:,:) = 0.e0 
     418      zqsro   (:,:) = 0.e0 
     419      zqsrics (:,:) = 0.e0 
     420      zqsrios (:,:) = 0.e0 
     421 
     422      DO jt = 1, jpintsr  
     423#   if defined key_vectopt_loop  
     424         DO ji = 1, jpij   
     425            zalbocsd(ji,1) = zalbocsd(ji,1) + zdlha   (ji,1) * zalbocs(ji,1,jt)   & 
     426               &                                             / MAX( 2.0 * zlsrise(ji,1) , zeps0 ) 
     427            zqsro   (ji,1) = zqsro   (ji,1) + zsqsro  (ji,1,jt) 
     428            zqsrics (ji,1) = zqsrics (ji,1) + zsqsrics(ji,1,jt) 
     429            zqsrios (ji,1) = zqsrios (ji,1) + zsqsrios(ji,1,jt) 
     430         END DO 
     431#  else 
     432         DO jj = 1, jpj 
     433            DO ji = 1, jpi   
     434               zalbocsd(ji,jj) = zalbocsd(ji,jj) + zdlha(ji,jj) * zalbocs(ji,jj,jt)   & 
     435                  &                                              / MAX( 2.0 * zlsrise(ji,jj) , zeps0 ) 
     436               zqsro  (ji,jj)  = zqsro   (ji,jj) + zsqsro  (ji,jj,jt) 
     437               zqsrics(ji,jj)  = zqsrics (ji,jj) + zsqsrics(ji,jj,jt) 
     438               zqsrios(ji,jj)  = zqsrios (ji,jj) + zsqsrios(ji,jj,jt) 
     439            END DO 
     440         END DO 
     441#  endif 
     442      END DO 
     443 
     444      DO jj = 1, jpj 
     445         DO ji = 1, jpi  
     446 
     447            !-------------------------------------------  
     448            !  Computation of shortwave radiation. 
     449            !------------------------------------------- 
     450 
     451            ! the solar heat flux absorbed at ocean and snow/ice surfaces 
     452            !------------------------------------------------------------ 
     453 
     454            ! For snow/ice  
     455            qsr_ice(ji,jj,jl) = ( zcatm1(ji,jj) * zqsrics(ji,jj) + catm(ji,jj) * zqsrios(ji,jj) ) / z2pi 
     456 
     457            ! Taking into account the ellipsity of the earth orbit 
     458            !----------------------------------------------------- 
     459 
     460            qsr_ice(ji,jj,jl) = qsr_ice(ji,jj,jl) * zdaycor 
     461            !--------------------------------------------------------------------------- 
     462            !   Computation of long-wave radiation  ( Berliand 1952 ; all latitudes ) 
     463            !--------------------------------------------------------------------------- 
     464 
     465            ! tempory variables 
     466            ztatm3         = ztatm(ji,jj) * ztatm(ji,jj) * ztatm(ji,jj) 
     467            ztatm4         = ztatm3 * ztatm(ji,jj) 
     468            z4tatm3        = 4. * ztatm3 
     469            zcldeff(ji,jj) = 1.0 - sbudyko(ji,jj) * catm(ji,jj) * catm(ji,jj)     
     470            ztaevbk        = ztatm4 * zcldeff(ji,jj) * ( 0.39 - 0.05 * zevsqr(ji,jj) )  
     471 
     472            !  Long-Wave for Ice 
     473            !---------------------- 
     474            zqlw_ice(ji,jj,jl) = - emic * stefan * ( ztaevbk + z4tatm3 * ( t_su(ji,jj,jl) - ztatm(ji,jj) ) )  
     475 
     476         END DO !ji 
     477      END DO !jj 
     478 
     479      END DO !jl 
     480 
     481      DO jj = 1, jpj 
     482         DO ji = 1, jpi  
     483 
     484            !  fraction of net shortwave radiation which is not absorbed in the  
     485            !  thin surface layer and penetrates inside the ice cover  
     486            !  ( Maykut and Untersteiner, 1971 ; Elbert anbd Curry, 1993 ) 
     487            !------------------------------------------------------------------ 
     488 
     489            fr1_i0(ji,jj) = 0.18  * zcatm1(ji,jj) + 0.35 * catm(ji,jj)  
     490            fr2_i0(ji,jj) = 0.82  * zcatm1(ji,jj) + 0.65 * catm(ji,jj) 
     491 
     492            ! the solar heat flux absorbed at ocean and snow/ice surfaces 
     493            !------------------------------------------------------------ 
     494            ! For ocean 
     495            qsr_oce(ji,jj) = srgamma * zcldcor(ji,jj) * zqsro(ji,jj) / z2pi 
     496            zinda          = SIGN( zone , -( -0.5 - zfrld(ji,jj) ) ) 
     497            zinda          = 1.0 - MAX( zzero , zinda ) 
     498            qsr_oce(ji,jj) = ( 1.- zinda ) * qsr_oce(ji,jj) 
     499 
     500            ! Taking into account the ellipsity of the earth orbit 
     501            !----------------------------------------------------- 
     502            qsr_oce(ji,jj) = qsr_oce(ji,jj) * zdaycor 
     503 
     504            !--------------------------------------------------------------------------- 
     505            !   Computation of long-wave radiation  ( Berliand 1952 ; all latitudes ) 
     506            !--------------------------------------------------------------------------- 
     507 
     508            ! tempory variables 
     509            ztatm3         = ztatm(ji,jj) * ztatm(ji,jj) * ztatm(ji,jj) 
     510            ztatm4         = ztatm3 * ztatm(ji,jj) 
     511            z4tatm3        = 4. * ztatm3 
     512            zcldeff(ji,jj) = 1.0 - sbudyko(ji,jj) * catm(ji,jj) * catm(ji,jj)     
     513            ztaevbk        = ztatm4 * zcldeff(ji,jj) * ( 0.39 - 0.05 * zevsqr(ji,jj) )  
     514 
     515            !  Long-Wave for Ocean 
     516            !----------------------- 
     517            zqlw_oce(ji,jj) = - emic * stefan * ( ztaevbk + z4tatm3 * ( psst  (ji,jj) - ztatm(ji,jj) ) )  
     518 
     519         END DO 
     520      END DO 
     521 
     522      !---------------------------------------- 
     523      !  Computation of turbulent heat fluxes  ( Latent and sensible )  
     524      !----------------------------------------         
     525      !CDIR NOVERRCHK 
     526      DO jj = 2 , jpjm1 
     527         !CDIR NOVERRCHK 
     528         DO  ji = 1, jpi 
     529 
     530            !  Turbulent heat fluxes over water 
     531            !---------------------------------- 
     532 
     533            ! zeso     : vapour pressure at saturation of ocean 
     534            ! zqsato   : humidity close to the ocean surface (at saturation) 
     535            zeso          =  611.0 * EXP ( 17.2693884 * ( psst(ji,jj) - rtt ) * tmask(ji,jj,1) / ( psst(ji,jj) - 35.86 ) ) 
     536            zqsato        = ( 0.622 * zeso ) / ( zpatm(ji,jj) - 0.378 * zeso ) 
     537 
     538            !  Drag coefficients from Large and Pond (1981,1982) 
     539            !--------------------------------------------------- 
     540     
     541            !  Stability parameters 
     542            zdteta         = psst(ji,jj) - ztatm(ji,jj) 
     543            zdeltaq        = zqatm(ji,jj) - zqsato 
     544            ztvmoy         = ztatm(ji,jj) * ( 1. + 2.2e-3 * ztatm(ji,jj) * zqatm(ji,jj) ) 
     545            zdenum         = MAX( vatm(ji,jj) * vatm(ji,jj) * ztvmoy, zeps ) 
     546            zdtetar        = zdteta / zdenum 
     547            ztvmoyr        = ztvmoy * ztvmoy * zdeltaq / zdenum 
     548             
     549            ! For stable atmospheric conditions 
     550            zobouks        = -70.0 * 10. * ( zdtetar + 3.2e-3 * ztvmoyr ) 
     551            zobouks        = MAX( zzero , zobouks ) 
     552            zpsims         = -7.0 * zobouks 
     553            zpsihs         =  zpsims 
     554            zpsils         =  zpsims 
     555 
     556            !  For unstable atmospheric conditions 
     557            zobouku        = -100.0 * 10.0 * ( zdtetar + 2.2e-3 * ztvmoyr ) 
     558            zobouku        = MIN( zzero , zobouku ) 
     559            zxins          = ( 1. - 16. * zobouku )**0.25 
     560            zlxins         = LOG( ( 1. + zxins * zxins ) / 2. ) 
     561            zpsimu         = 2. * LOG( ( 1 + zxins ) / 2. )  + zlxins - 2. * ATAN( zxins ) + zpis2 
     562            zpsihu         = 2. * zlxins 
     563            zpsilu         = zpsihu 
     564 
     565            ! computation of intermediate values 
     566            zstab          = MAX( zzero , SIGN( zone , zdteta ) ) 
     567            zpsim          = zstab * zpsimu + (1.0 - zstab ) * zpsims 
     568            zpsih          = zstab * zpsihu + (1.0 - zstab ) * zpsihs 
     569            zpsil          = zpsih 
     570             
     571            zvatmg         = MAX( 0.032 * 1.5e-3 * vatm(ji,jj) * vatm(ji,jj) / grav, zeps ) 
     572 
     573            zcmn           = vkarmn / LOG ( 10. / zvatmg ) 
     574            zcmn2          = zcmn * zcmn 
     575            zchn           = 0.0327 * zcmn 
     576            zcln           = 0.0346 * zcmn 
     577            zcmcmn         = 1 / ( 1 - zcmn * zpsim / vkarmn ) 
     578            zchcm          = zcmcmn / ( 1 - zchn * zpsih / ( vkarmn * zcmn ) ) 
     579            zclcm          = zchcm 
     580 
     581 
     582            !  Transfer cofficient zcsho and zcleo over ocean according to Large and Pond (1981,1982) 
     583            !--------------------------------------------------------------  
     584            zcsho          = zchn * zchcm 
     585            zcleo          = zcln * zclcm  
     586 
     587 
     588            !   Computation of sensible and latent fluxes over Ocean  
     589            !---------------------------------------------------------------- 
     590 
     591            !  computation of intermediate values 
     592            zrhova         = zrhoa(ji,jj) * vatm(ji,jj) 
     593            zrhovacsho     = zrhova * zcsho 
     594            zrhovacleo     = zrhova * zcleo 
     595 
     596            ! sensible heat flux 
     597            zqsb_oce(ji,jj) = zrhovacsho * 1004.0  * ( psst(ji,jj) - ztatm(ji,jj) )   
     598          
     599            !  latent heat flux  
     600            zqla_oce(ji,jj) = MAX(0.e0, zrhovacleo * 2.5e+06 * ( zqsato      - zqatm(ji,jj) ) ) 
     601                
     602            !  Calculate evaporation over water. (kg/m2/s) 
     603            !------------------------------------------------- 
     604            evap(ji,jj)    = zqla_oce(ji,jj) / cevap 
     605                
     606         END DO !ji 
     607      END DO !jj 
     608 
     609      DO jl = 1, jpl 
     610      !CDIR NOVERRCHK 
     611      DO jj = 2 , jpjm1 
     612         !CDIR NOVERRCHK 
     613         DO ji = 1, jpi 
     614                
     615            !  Turbulent heat fluxes over snow/ice. 
     616            !-------------------------------------------------- 
     617             
     618            !  zesi     : vapour pressure at saturation of ice 
     619            !  zqsati   : humidity close to the ice surface (at saturation) 
     620            zesi           =  611.0 * EXP ( 21.8745587 * tmask(ji,jj,1)   &   ! tmask needed to avoid overflow in the exponential 
     621               &                                       * ( t_su(ji,jj,jl) - rtt )/ ( t_su(ji,jj,jl)- 7.66 ) ) 
     622            zqsati         = ( 0.622 * zesi ) / ( zpatm(ji,jj) - 0.378 * zesi ) 
     623                
     624            !  computation of intermediate values 
     625            zticemb        = t_su(ji,jj,jl) - 7.66 
     626            zticemb2       = zticemb * zticemb   
     627            ztice3         = t_su(ji,jj,jl) * t_su(ji,jj,jl) * t_su(ji,jj,jl) 
     628            zqsati2        = zqsati * zqsati 
     629            zesi2          = zesi * zesi 
     630            zdesidt        = zesi * ( 9.5 * LOG( 10.0 ) * ( rtt - 7.66 )  / zticemb2 ) 
     631                
     632            !  Transfer cofficient zcshi and zclei over ice. Assumed to be constant Parkinson 1979 ; Maykut 1982 
     633            !-------------------------------------------------------------------- 
     634            zcshi          = 1.75e-03 
     635            zclei          = zcshi 
     636                
     637            !  Computation of sensible and latent fluxes over ice 
     638            !---------------------------------------------------------------- 
     639                
     640            !  computation of intermediate values 
     641            zrhova          = zrhoa(ji,jj) * vatm(ji,jj) 
     642            zrhovaclei      = zrhova * zcshi * 2.834e+06 
     643            zrhovacshi      = zrhova * zclei * 1004.0 
     644             
     645            !  sensible heat flux 
     646            zqsb_ice(ji,jj,jl) = zrhovacshi * ( t_su(ji,jj,jl) - ztatm(ji,jj) ) 
     647             
     648            !  latent heat flux  
     649            zqla_ice(ji,jj,jl) = zrhovaclei * ( zqsati        - zqatm(ji,jj) ) 
     650            qla_ice (ji,jj,jl) = MAX(0.e0, zqla_ice(ji,jj,jl) ) 
     651               
     652            !  Computation of sensitivity of non solar fluxes (dQ/dT) 
     653            !--------------------------------------------------------------- 
     654                
     655            !  computation of long-wave, sensible and latent flux sensitivity 
     656            zdqlw_ice       = 4.0 * emic * stefan * ztice3 
     657            zdqsb_ice       = zrhovacshi 
     658            zdqla_ice       = zrhovaclei * ( zdesidt * ( zqsati2 / zesi2 ) * ( zpatm(ji,jj) / 0.622 ) )    
     659             
     660            !  total non solar sensitivity 
     661            dqns_ice(ji,jj,jl) = -( zdqlw_ice + zdqsb_ice + zdqla_ice )  
     662             
     663            ! latent flux sensitivity 
     664            dqla_ice(ji,jj,jl) = zdqla_ice 
     665             
     666         END DO 
     667      END DO 
     668      END DO !jl 
     669 
     670      ! total non solar heat flux over ice 
     671      qnsr_ice(:,:,:) = zqlw_ice(:,:,:) - zqsb_ice(:,:,:) - zqla_ice(:,:,:) 
     672      ! total non solar heat flux over water  
     673      qnsr_oce(:,:) = zqlw_oce(:,:) - zqsb_oce(:,:) - zqla_oce(:,:) 
     674 
     675      ! solid precipitations ( kg/m2/day -> kg/m2/s) 
     676      tprecip(:,:) = tprecip  (:,:) / rday  
     677      ! snow  precipitations ( kg/m2/day -> kg/m2/s) 
     678      sprecip(:,:) = sprecip  (:,:) / rday   
     679 
     680      CALL lbc_lnk( qsr_oce (:,:) , 'T', 1. ) 
     681      CALL lbc_lnk( qnsr_oce(:,:) , 'T', 1. ) 
     682      CALL lbc_lnk( fr1_i0  (:,:) , 'T', 1. ) 
     683      CALL lbc_lnk( fr2_i0  (:,:) , 'T', 1. ) 
     684      CALL lbc_lnk( tprecip (:,:) , 'T', 1. ) 
     685      CALL lbc_lnk( sprecip (:,:) , 'T', 1. ) 
     686      CALL lbc_lnk( evap    (:,:) , 'T', 1. ) 
     687      DO jl = 1, jpl 
     688         CALL lbc_lnk( qsr_ice (:,:,jl) , 'T', 1. ) 
     689         CALL lbc_lnk( qnsr_ice(:,:,jl) , 'T', 1. ) 
     690         CALL lbc_lnk( dqns_ice(:,:,jl) , 'T', 1. ) 
     691         CALL lbc_lnk( qla_ice (:,:,jl) , 'T', 1. ) 
     692         CALL lbc_lnk( dqla_ice(:,:,jl) , 'T', 1. ) 
     693      END DO 
     694 
     695      qsr_oce (:,:) = qsr_oce (:,:)*tmask(:,:,1) 
     696      qnsr_oce(:,:) = qnsr_oce(:,:)*tmask(:,:,1) 
     697      DO jl = 1, jpl 
     698         qsr_ice (:,:,jl) = qsr_ice (:,:,jl)*tmask(:,:,1) 
     699         qnsr_ice(:,:,jl) = qnsr_ice(:,:,jl)*tmask(:,:,1) 
     700         qla_ice (:,:,jl) = qla_ice (:,:,jl)*tmask(:,:,1) 
     701         dqns_ice(:,:,jl) = dqns_ice(:,:,jl)*tmask(:,:,1) 
     702         dqla_ice(:,:,jl) = dqla_ice(:,:,jl)*tmask(:,:,1) 
     703      END DO 
     704      fr1_i0  (:,:) = fr1_i0  (:,:)*tmask(:,:,1) 
     705      fr2_i0  (:,:) = fr2_i0  (:,:)*tmask(:,:,1) 
     706      tprecip (:,:) = tprecip (:,:)*tmask(:,:,1) 
     707      sprecip (:,:) = sprecip (:,:)*tmask(:,:,1) 
     708      evap    (:,:) = evap    (:,:)*tmask(:,:,1) 
     709 
     710 
     711   END SUBROUTINE flx_blk 
     712 
     713 
     714#elif defined key_lim2 
    81715 
    82716   SUBROUTINE flx_blk( psst ) 
     
    212846      !!--------------------------------------------------------------------- 
    213847 
    214       !--------------------- 
    215       !   Initilization    ! 
    216       !--------------------- 
    217 #if ! defined key_lim3 
    218       tn_ice(:,:) = psst(:,:) 
    219 #endif 
    220  
    221848      !  Determine cloud optical depths as a function of latitude (Chou et al., 1981). 
    222849      !  and the correction factor for taking into account  the effect of clouds  
     
    6801307 
    6811308   END SUBROUTINE flx_blk 
     1309#endif 
    6821310 
    6831311 
  • branches/dev_002_LIM/NEMO/OPA_SRC/SBC/ocesbc.F90

    r826 r827  
    6666#if defined key_lim3 || defined key_lim2 
    6767   !!---------------------------------------------------------------------- 
    68    !!   'key_lim3'     :                              LIM2 sea-ice model 
    69    !!   'key_lim2'        :                              LIM3 sea-ice model 
     68   !!   'key_lim3'     :                                 LIM2 sea-ice model 
     69   !!   'key_lim2'     :                                 LIM3 sea-ice model 
    7070   !!---------------------------------------------------------------------- 
    7171# if defined key_coupled 
     
    155155            DO ji = 1, fs_jpim1   ! vertor opt. 
    156156               ztx   = 0.5 * ( freeze(ji+1,jj) + freeze(ji+1,jj+1) ) 
     157               zty   = 0.5 * ( freeze(ji,jj+1) + freeze(ji+1,jj+1) ) 
     158#if defined key_lim3 
     159               ztaux = ftaux(ji,jj) 
     160               ztauy = ftauy(ji,jj) 
     161#elif defined key_lim2 
    157162               ztaux = 0.5 * ( ftaux (ji+1,jj) + ftaux (ji+1,jj+1) ) 
     163               ztauy = 0.5 * ( ftauy (ji,jj+1) + ftauy (ji+1,jj+1) ) 
     164#endif 
    158165               taux(ji,jj) = (1.0-ztx) * taux(ji,jj) + ztx * ztaux 
    159  
    160                zty   = 0.5 * ( freeze(ji,jj+1) + freeze(ji+1,jj+1) ) 
    161                ztauy = 0.5 * ( ftauy (ji,jj+1) + ftauy (ji+1,jj+1) ) 
    162166               tauy(ji,jj) = (1.0-zty) * tauy(ji,jj) + zty * ztauy 
    163167            END DO 
     
    260264               ztx         = MAX( freezn(ji,jj), freezn(ji,jj+1) )   ! ice/ocean indicator at U- and V-points 
    261265               zty         = MAX( freezn(ji,jj), freezn(ji+1,jj) ) 
    262                ztaux       = 0.5 *( ftaux(ji+1,jj) + ftaux(ji+1,jj+1) ) ! ice-ocean stress at U- and V-points 
    263                ztauy       = 0.5 *( ftauy(ji,jj+1) + ftauy(ji+1,jj+1) ) 
     266#if defined key_lim3 
     267               ztaux = ftaux(ji,jj) 
     268               ztauy = ftauy(ji,jj) 
     269#elif defined key_lim2 
     270               ztaux = 0.5 * ( ftaux (ji+1,jj) + ftaux (ji+1,jj+1) ) 
     271               ztauy = 0.5 * ( ftauy (ji,jj+1) + ftauy (ji+1,jj+1) ) 
     272#endif 
    264273               taux(ji,jj) = (1.-ztx) * taux(ji,jj) + ztx * ztaux    ! stress at the ocean surface 
    265274               tauy(ji,jj) = (1.-zty) * tauy(ji,jj) + zty * ztauy 
  • branches/dev_002_LIM/NEMO/OPA_SRC/ice_oce.F90

    r826 r827  
    2525# if defined  key_lim2 
    2626   LOGICAL, PUBLIC, PARAMETER ::   lk_lim2        = .TRUE.    !: LIM2 ice model 
    27    LOGICAL, PUBLIC, PARAMETER ::   lk_ice_lim     = .FALSE.   !: LIM3 ice model 
     27   LOGICAL, PUBLIC, PARAMETER ::   lk_lim3        = .FALSE.   !: LIM3 ice model 
    2828# else 
    2929   LOGICAL, PUBLIC, PARAMETER ::   lk_lim2        = .FALSE.   !: LIM2 ice model 
    30    LOGICAL, PUBLIC, PARAMETER ::   lk_ice_lim     = .TRUE.    !: LIM3 ice model 
     30   LOGICAL, PUBLIC, PARAMETER ::   lk_lim3        = .TRUE.    !: LIM3 ice model 
    3131# endif 
    3232 
     
    3737   REAL(wp), PUBLIC, DIMENSION(jpiglo,jpjglo) ::   &  !: cumulated fields 
    3838      fqsr_oce ,      &   !: Net short wave heat flux on free ocean  
    39       fqsr_ice ,      &   !: Net short wave het flux on sea ice  
     39      fqsr_ice ,      &   !: Net short wave heat flux on sea ice  
    4040      fqnsr_oce,      &   !: Net longwave heat flux on free ocean 
    4141      fqnsr_ice,      &   !: Net longwave heat flux on sea ice 
     
    5454      ftaux , ftauy  , &  !: wind stresses 
    5555      gtaux , gtauy       !: wind stresses 
     56 
     57# if defined key_lim3 
     58   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &  !: field exchanges with ice model to ocean 
     59      catm_ice       , &  !: cloud cover 
     60      tatm_ice       , &  !: air temperature 
     61      icethi              !: icethickness 
     62# endif 
    5663    
    5764   REAL(wp), PUBLIC ::   &  !: 
     
    6471   !!---------------------------------------------------------------------- 
    6572   LOGICAL, PUBLIC, PARAMETER ::   lk_lim2        = .FALSE.  !: No LIM 2.0 ice model 
    66    LOGICAL, PUBLIC, PARAMETER ::   lk_ice_lim     = .FALSE.  !: No LIM 3.0 ice model 
     73   LOGICAL, PUBLIC, PARAMETER ::   lk_lim3        = .FALSE.  !: No LIM 3.0 ice model 
    6774#endif 
    6875 
  • branches/dev_002_LIM/NEMO/OPA_SRC/phycst.F90

    r719 r827  
    4646      rtt      = 273.16_wp  ,  &  !: triple point of temperature (Kelvin) 
    4747      rt0      = 273.15_wp  ,  &  !: freezing point of water (Kelvin) 
     48#if defined key_lim3 
     49      rt0_snow = 273.16_wp  ,  &  !: melting point of snow  (Kelvin) 
     50      rt0_ice  = 273.16_wp  ,  &  !: melting point of ice   (Kelvin) 
     51#elif defined key_lim2 
    4852      rt0_snow = 273.15_wp  ,  &  !: melting point of snow  (Kelvin) 
    4953      rt0_ice  = 273.05_wp  ,  &  !: melting point of ice   (Kelvin) 
     54#endif 
    5055      rau0     = 1020._wp   ,  &  !: volumic mass of reference (kg/m3) 
    5156      rauw     = 1000._wp   ,  &  !: density of pure water (kg/m3) 
     
    5459 
    5560   REAL(wp), PUBLIC ::            &  !: 
     61#if defined key_lim3 
     62      rcdsn   =   0.31_wp     ,   &  !: thermal conductivity of snow 
     63      rcdic   =   2.034396_wp ,   &  !: thermal conductivity of fresh ice 
     64      cpic    = 2067.0        ,   & 
     65      ! add the following lines 
     66      lsub    = 2.834e+6      ,   &  !: pure ice latent heat of sublimation (J.kg-1) 
     67      lfus    = 0.334e+6      ,   &  !: latent heat of fusion of fresh ice   (J.kg-1) 
     68      rhoic   = 917._wp       ,   &  !: density of sea ice (kg/m3) 
     69      tmut    =   0.054       ,   &  !: decrease of seawater meltpoint with salinity 
     70#elif defined key_lim2 
    5671      rcdsn   =   0.22_wp     ,   &  !: conductivity of the snow 
    5772      rcdic   =   2.034396_wp ,   &  !: conductivity of the ice 
     
    6277      xsn     =   2.8e+6      ,   &  !: latent heat of sublimation of snow 
    6378      rhoic   = 900._wp       ,   &  !: density of sea ice (kg/m3) 
     79#endif 
    6480      rhosn   = 330._wp       ,   &  !: density of snow (kg/m3) 
    6581      emic    =   0.97_wp     ,   &  !: emissivity of snow or ice 
     
    169185         WRITE(numout,*) '          thermal conductivity of the snow          = ', rcdsn   , ' J/s/m/K' 
    170186         WRITE(numout,*) '          thermal conductivity of the ice           = ', rcdic   , ' J/s/m/K' 
     187#if defined key_lim3 
     188         WRITE(numout,*) '          fresh ice specific heat                   = ', cpic    , ' J/kg/K' 
     189         WRITE(numout,*) '          latent heat of fusion of fresh ice / snow = ', lfus    , ' J/kg' 
     190         WRITE(numout,*) '          latent heat of subl.  of fresh ice / snow = ', lsub    , ' J/kg' 
     191#elif defined key_lim2 
    171192         WRITE(numout,*) '          density times specific heat for snow      = ', rcpsn   , ' J/m^3/K'  
    172193         WRITE(numout,*) '          density times specific heat for ice       = ', rcpic   , ' J/m^3/K' 
     
    174195         WRITE(numout,*) '          volumetric latent heat fusion of snow     = ', xlsn    , ' J/m'  
    175196         WRITE(numout,*) '          latent heat of sublimation of snow        = ', xsn     , ' J/kg'  
     197#endif 
    176198         WRITE(numout,*) '          density of sea ice                        = ', rhoic   , ' kg/m^3' 
    177199         WRITE(numout,*) '          density of snow                           = ', rhosn   , ' kg/m^3' 
Note: See TracChangeset for help on using the changeset viewer.