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 14037 for NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/SBC/sbccpl.F90 – NEMO

Ignore:
Timestamp:
2020-12-03T12:20:38+01:00 (3 years ago)
Author:
ayoung
Message:

Updated to trunk at 14020. Sette tests passed with change of results for configurations with non-linear ssh. Ticket #2506.

Location:
NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG

    • Property svn:externals
      •  

        old new  
        88 
        99# SETTE 
        10 ^/utils/CI/sette@13292        sette 
         10^/utils/CI/sette_wave@13990         sette 
  • NEMO/branches/2020/dev_r13333_KERNEL-08_techene_gm_HPG_SPG/src/OCE/SBC/sbccpl.F90

    r13295 r14037  
    88   !!            3.1  ! 2009_02  (G. Madec, S. Masson, E. Maisonave, A. Caubel) generic coupled interface 
    99   !!            3.4  ! 2011_11  (C. Harris) more flexibility + multi-category fields 
     10   !!            4.2  ! 2020-12  (G. Madec, E. Clementi)  wave coupling updates 
    1011   !!---------------------------------------------------------------------- 
    1112 
     
    4142#endif 
    4243#if defined key_si3 
    43    USE icethd_dh      ! for CALL ice_thd_snwblow 
     44   USE icevar         ! for CALL ice_var_snwblow 
    4445#endif 
    4546   ! 
     
    4849   USE lib_mpp        ! distribued memory computing library 
    4950   USE lbclnk         ! ocean lateral boundary conditions (or mpp link) 
     51 
     52#if defined key_oasis3  
     53   USE mod_oasis, ONLY : OASIS_Sent, OASIS_ToRest, OASIS_SentOut, OASIS_ToRestOut  
     54#endif  
    5055 
    5156   IMPLICIT NONE 
     
    102107   INTEGER, PARAMETER ::   jpr_fraqsr = 42   ! fraction of solar net radiation absorbed in the first ocean level 
    103108   INTEGER, PARAMETER ::   jpr_mslp   = 43   ! mean sea level pressure  
    104    INTEGER, PARAMETER ::   jpr_hsig   = 44   ! Hsig  
    105    INTEGER, PARAMETER ::   jpr_phioc  = 45   ! Wave=>ocean energy flux  
    106    INTEGER, PARAMETER ::   jpr_sdrftx = 46   ! Stokes drift on grid 1  
    107    INTEGER, PARAMETER ::   jpr_sdrfty = 47   ! Stokes drift on grid 2  
     109   !**  surface wave coupling  ** 
     110   INTEGER, PARAMETER ::   jpr_hsig   = 44   ! Hsig 
     111   INTEGER, PARAMETER ::   jpr_phioc  = 45   ! Wave=>ocean energy flux 
     112   INTEGER, PARAMETER ::   jpr_sdrftx = 46   ! Stokes drift on grid 1 
     113   INTEGER, PARAMETER ::   jpr_sdrfty = 47   ! Stokes drift on grid 2 
    108114   INTEGER, PARAMETER ::   jpr_wper   = 48   ! Mean wave period 
    109115   INTEGER, PARAMETER ::   jpr_wnum   = 49   ! Mean wavenumber 
    110    INTEGER, PARAMETER ::   jpr_tauwoc = 50   ! Stress fraction adsorbed by waves 
     116   INTEGER, PARAMETER ::   jpr_wstrf = 50   ! Stress fraction adsorbed by waves 
    111117   INTEGER, PARAMETER ::   jpr_wdrag  = 51   ! Neutral surface drag coefficient 
    112    INTEGER, PARAMETER ::   jpr_isf    = 52 
    113    INTEGER, PARAMETER ::   jpr_icb    = 53 
    114    INTEGER, PARAMETER ::   jpr_wfreq  = 54   ! Wave peak frequency 
    115    INTEGER, PARAMETER ::   jpr_tauwx  = 55   ! x component of the ocean stress from waves 
    116    INTEGER, PARAMETER ::   jpr_tauwy  = 56   ! y component of the ocean stress from waves 
    117    INTEGER, PARAMETER ::   jpr_ts_ice = 57   ! Sea ice surface temp 
    118  
    119    INTEGER, PARAMETER ::   jprcv      = 57   ! total number of fields received   
     118   INTEGER, PARAMETER ::   jpr_charn  = 52   ! Chranock coefficient 
     119   INTEGER, PARAMETER ::   jpr_twox   = 53   ! wave to ocean momentum flux 
     120   INTEGER, PARAMETER ::   jpr_twoy   = 54   ! wave to ocean momentum flux 
     121   INTEGER, PARAMETER ::   jpr_tawx   = 55   ! net wave-supported stress 
     122   INTEGER, PARAMETER ::   jpr_tawy   = 56   ! net wave-supported stress 
     123   INTEGER, PARAMETER ::   jpr_bhd    = 57   ! Bernoulli head. waves' induced surface pressure 
     124   INTEGER, PARAMETER ::   jpr_tusd   = 58   ! zonal stokes transport 
     125   INTEGER, PARAMETER ::   jpr_tvsd   = 59   ! meridional stokes tranmport 
     126   INTEGER, PARAMETER ::   jpr_isf    = 60 
     127   INTEGER, PARAMETER ::   jpr_icb    = 61 
     128   INTEGER, PARAMETER ::   jpr_ts_ice = 62   ! Sea ice surface temp 
     129 
     130   INTEGER, PARAMETER ::   jprcv      = 62   ! total number of fields received   
    120131 
    121132   INTEGER, PARAMETER ::   jps_fice   =  1   ! ice fraction sent to the atmosphere 
     
    152163   INTEGER, PARAMETER ::   jps_wlev   = 32   ! water level  
    153164   INTEGER, PARAMETER ::   jps_fice1  = 33   ! first-order ice concentration (for semi-implicit coupling of atmos-ice fluxes) 
    154    INTEGER, PARAMETER ::   jps_a_p    = 34   ! meltpond area 
     165   INTEGER, PARAMETER ::   jps_a_p    = 34   ! meltpond area fraction 
    155166   INTEGER, PARAMETER ::   jps_ht_p   = 35   ! meltpond thickness 
    156167   INTEGER, PARAMETER ::   jps_kice   = 36   ! sea ice effective conductivity 
     
    159170 
    160171   INTEGER, PARAMETER ::   jpsnd      = 38   ! total number of fields sent  
     172 
     173#if ! defined key_oasis3  
     174   ! Dummy variables to enable compilation when oasis3 is not being used  
     175   INTEGER                    ::   OASIS_Sent        = -1  
     176   INTEGER                    ::   OASIS_SentOut     = -1  
     177   INTEGER                    ::   OASIS_ToRest      = -1  
     178   INTEGER                    ::   OASIS_ToRestOut   = -1  
     179#endif  
    161180 
    162181   !                                  !!** namelist namsbc_cpl ** 
     
    172191      &             sn_snd_thick1, sn_snd_cond, sn_snd_mpnd , sn_snd_sstfrz, sn_snd_ttilyr 
    173192   !                                   ! Received from the atmosphere 
    174    TYPE(FLD_C) ::   sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_tauw, sn_rcv_dqnsdt, sn_rcv_qsr,  & 
     193   TYPE(FLD_C) ::   sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr,  & 
    175194      &             sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf, sn_rcv_ts_ice 
    176195   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_mslp, sn_rcv_icb, sn_rcv_isf 
    177    ! Send to waves  
     196   !                                   ! Send to waves  
    178197   TYPE(FLD_C) ::   sn_snd_ifrac, sn_snd_crtw, sn_snd_wlev  
    179    ! Received from waves  
    180    TYPE(FLD_C) ::   sn_rcv_hsig, sn_rcv_phioc, sn_rcv_sdrfx, sn_rcv_sdrfy, sn_rcv_wper, sn_rcv_wnum, sn_rcv_tauwoc, & 
    181                     sn_rcv_wdrag, sn_rcv_wfreq 
     198   !                                   ! Received from waves  
     199   TYPE(FLD_C) ::   sn_rcv_hsig, sn_rcv_phioc, sn_rcv_sdrfx, sn_rcv_sdrfy, sn_rcv_wper, sn_rcv_wnum, & 
     200      &             sn_rcv_wstrf, sn_rcv_wdrag, sn_rcv_charn, sn_rcv_taw, sn_rcv_bhd, sn_rcv_tusd, sn_rcv_tvsd 
    182201   !                                   ! Other namelist parameters 
    183202   INTEGER     ::   nn_cplmodel           ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 
    184203   LOGICAL     ::   ln_usecplmask         !  use a coupling mask file to merge data received from several models 
    185204                                          !   -> file cplmask.nc with the float variable called cplmask (jpi,jpj,nn_cplmodel) 
     205   LOGICAL     ::   ln_scale_ice_flux     !  use ice fluxes that are already "ice weighted" ( i.e. multiplied ice concentration)  
     206 
    186207   TYPE ::   DYNARR      
    187208      REAL(wp), POINTER, DIMENSION(:,:,:) ::   z3    
     
    191212 
    192213   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:) ::   alb_oce_mix    ! ocean albedo sent to atmosphere (mix clear/overcast sky) 
     214#if defined key_si3 || defined key_cice 
     215   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   a_i_last_couple !: Ice fractional area at last coupling time 
     216#endif 
    193217 
    194218   REAL(wp) ::   rpref = 101000._wp   ! reference atmospheric pressure[N/m2]  
     
    211235      !!             ***  FUNCTION sbc_cpl_alloc  *** 
    212236      !!---------------------------------------------------------------------- 
    213       INTEGER :: ierr(4) 
     237      INTEGER :: ierr(5) 
    214238      !!---------------------------------------------------------------------- 
    215239      ierr(:) = 0 
     
    221245#endif 
    222246      ALLOCATE( xcplmask(jpi,jpj,0:nn_cplmodel) , STAT=ierr(3) ) 
    223       ! 
    224       IF( .NOT. ln_apr_dyn ) ALLOCATE( ssh_ib(jpi,jpj), ssh_ibb(jpi,jpj), apr(jpi, jpj), STAT=ierr(4) )  
     247#if defined key_si3 || defined key_cice 
     248      ALLOCATE( a_i_last_couple(jpi,jpj,jpl) , STAT=ierr(4) ) 
     249#endif 
     250      ! 
     251      IF( .NOT. ln_apr_dyn ) ALLOCATE( ssh_ib(jpi,jpj), ssh_ibb(jpi,jpj), apr(jpi, jpj), STAT=ierr(5) ) 
    225252 
    226253      sbc_cpl_alloc = MAXVAL( ierr ) 
     
    249276      REAL(wp), DIMENSION(jpi,jpj) ::   zacs, zaos 
    250277      !! 
    251       NAMELIST/namsbc_cpl/  sn_snd_temp  , sn_snd_alb   , sn_snd_thick, sn_snd_crt   , sn_snd_co2  ,   &  
     278      NAMELIST/namsbc_cpl/  nn_cplmodel  , ln_usecplmask, nn_cats_cpl , ln_scale_ice_flux,             & 
     279         &                  sn_snd_temp  , sn_snd_alb   , sn_snd_thick, sn_snd_crt   , sn_snd_co2   ,  &  
    252280         &                  sn_snd_ttilyr, sn_snd_cond  , sn_snd_mpnd , sn_snd_sstfrz, sn_snd_thick1,  &  
    253          &                  sn_snd_ifrac , sn_snd_crtw  , sn_snd_wlev , sn_rcv_hsig  , sn_rcv_phioc,   &  
    254          &                  sn_rcv_w10m  , sn_rcv_taumod, sn_rcv_tau  , sn_rcv_dqnsdt, sn_rcv_qsr  ,   &  
    255          &                  sn_rcv_sdrfx , sn_rcv_sdrfy , sn_rcv_wper , sn_rcv_wnum  , sn_rcv_tauwoc,  & 
     281         &                  sn_snd_ifrac , sn_snd_crtw  , sn_snd_wlev , sn_rcv_hsig  , sn_rcv_phioc ,  &  
     282         &                  sn_rcv_w10m  , sn_rcv_taumod, sn_rcv_tau  , sn_rcv_dqnsdt, sn_rcv_qsr   ,  &  
     283         &                  sn_rcv_sdrfx , sn_rcv_sdrfy , sn_rcv_wper , sn_rcv_wnum  , sn_rcv_wstrf ,  & 
     284         &                  sn_rcv_charn , sn_rcv_taw   , sn_rcv_bhd  , sn_rcv_tusd  , sn_rcv_tvsd,    & 
    256285         &                  sn_rcv_wdrag , sn_rcv_qns   , sn_rcv_emp  , sn_rcv_rnf   , sn_rcv_cal  ,   & 
    257          &                  sn_rcv_iceflx, sn_rcv_co2   , nn_cplmodel , ln_usecplmask, sn_rcv_mslp ,   & 
    258          &                  sn_rcv_icb   , sn_rcv_isf   , sn_rcv_wfreq , sn_rcv_tauw, nn_cats_cpl  ,   & 
    259          &                  sn_rcv_ts_ice 
     286         &                  sn_rcv_iceflx, sn_rcv_co2   , sn_rcv_icb  , sn_rcv_isf   , sn_rcv_ts_ice  
    260287 
    261288      !!--------------------------------------------------------------------- 
     
    278305      ENDIF 
    279306      IF( lwp .AND. ln_cpl ) THEN                        ! control print 
     307         WRITE(numout,*)'  nn_cplmodel                         = ', nn_cplmodel 
     308         WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask 
     309         WRITE(numout,*)'  ln_scale_ice_flux                   = ', ln_scale_ice_flux 
     310         WRITE(numout,*)'  nn_cats_cpl                         = ', nn_cats_cpl 
    280311         WRITE(numout,*)'  received fields (mutiple ice categogies)' 
    281312         WRITE(numout,*)'      10m wind module                 = ', TRIM(sn_rcv_w10m%cldes  ), ' (', TRIM(sn_rcv_w10m%clcat  ), ')' 
     
    295326         WRITE(numout,*)'      sea ice heat fluxes             = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')' 
    296327         WRITE(numout,*)'      atm co2                         = ', TRIM(sn_rcv_co2%cldes   ), ' (', TRIM(sn_rcv_co2%clcat   ), ')' 
     328         WRITE(numout,*)'      Sea ice surface skin temperature= ', TRIM(sn_rcv_ts_ice%cldes), ' (', TRIM(sn_rcv_ts_ice%clcat), ')' 
     329         WRITE(numout,*)'      surface waves:' 
    297330         WRITE(numout,*)'      significant wave heigth         = ', TRIM(sn_rcv_hsig%cldes  ), ' (', TRIM(sn_rcv_hsig%clcat  ), ')'  
    298331         WRITE(numout,*)'      wave to oce energy flux         = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')'  
     
    301334         WRITE(numout,*)'      Mean wave period                = ', TRIM(sn_rcv_wper%cldes  ), ' (', TRIM(sn_rcv_wper%clcat  ), ')'  
    302335         WRITE(numout,*)'      Mean wave number                = ', TRIM(sn_rcv_wnum%cldes  ), ' (', TRIM(sn_rcv_wnum%clcat  ), ')'  
    303          WRITE(numout,*)'      Wave peak frequency             = ', TRIM(sn_rcv_wfreq%cldes ), ' (', TRIM(sn_rcv_wfreq%clcat ), ')' 
    304          WRITE(numout,*)'      Stress frac adsorbed by waves   = ', TRIM(sn_rcv_tauwoc%cldes), ' (', TRIM(sn_rcv_tauwoc%clcat ), ')'  
    305          WRITE(numout,*)'      Stress components by waves      = ', TRIM(sn_rcv_tauw%cldes  ), ' (', TRIM(sn_rcv_tauw%clcat  ), ')' 
     336         WRITE(numout,*)'      Stress frac adsorbed by waves   = ', TRIM(sn_rcv_wstrf%cldes ), ' (', TRIM(sn_rcv_wstrf%clcat ), ')' 
    306337         WRITE(numout,*)'      Neutral surf drag coefficient   = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')'  
    307          WRITE(numout,*)'      Sea ice surface skin temperature= ', TRIM(sn_rcv_ts_ice%cldes), ' (', TRIM(sn_rcv_ts_ice%clcat), ')'  
     338         WRITE(numout,*)'      Charnock coefficient            = ', TRIM(sn_rcv_charn%cldes ), ' (', TRIM(sn_rcv_charn%clcat ), ')' 
    308339         WRITE(numout,*)'  sent fields (multiple ice categories)' 
    309340         WRITE(numout,*)'      surface temperature             = ', TRIM(sn_snd_temp%cldes  ), ' (', TRIM(sn_snd_temp%clcat  ), ')' 
     
    326357         WRITE(numout,*)'                      - orientation   = ', sn_snd_crtw%clvor  
    327358         WRITE(numout,*)'                      - mesh          = ', sn_snd_crtw%clvgrd  
    328          WRITE(numout,*)'  nn_cplmodel                         = ', nn_cplmodel 
    329          WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask 
    330          WRITE(numout,*)'  nn_cats_cpl                         = ', nn_cats_cpl 
    331       ENDIF 
    332  
     359      ENDIF 
     360      IF( lwp .AND. ln_wave) THEN                        ! control print 
     361      WRITE(numout,*)'      surface waves:' 
     362         WRITE(numout,*)'      Significant wave heigth         = ', TRIM(sn_rcv_hsig%cldes  ), ' (', TRIM(sn_rcv_hsig%clcat  ), ')' 
     363         WRITE(numout,*)'      Wave to oce energy flux         = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')' 
     364         WRITE(numout,*)'      Surface Stokes drift grid u     = ', TRIM(sn_rcv_sdrfx%cldes ), ' (', TRIM(sn_rcv_sdrfx%clcat ), ')' 
     365         WRITE(numout,*)'      Surface Stokes drift grid v     = ', TRIM(sn_rcv_sdrfy%cldes ), ' (', TRIM(sn_rcv_sdrfy%clcat ), ')' 
     366         WRITE(numout,*)'      Mean wave period                = ', TRIM(sn_rcv_wper%cldes  ), ' (', TRIM(sn_rcv_wper%clcat  ), ')' 
     367         WRITE(numout,*)'      Mean wave number                = ', TRIM(sn_rcv_wnum%cldes  ), ' (', TRIM(sn_rcv_wnum%clcat  ), ')' 
     368         WRITE(numout,*)'      Stress frac adsorbed by waves   = ', TRIM(sn_rcv_wstrf%cldes ), ' (', TRIM(sn_rcv_wstrf%clcat ), ')' 
     369         WRITE(numout,*)'      Neutral surf drag coefficient   = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')' 
     370         WRITE(numout,*)'      Charnock coefficient            = ', TRIM(sn_rcv_charn%cldes ), ' (', TRIM(sn_rcv_charn%clcat ), ')' 
     371         WRITE(numout,*)' Transport associated to Stokes drift grid u = ', TRIM(sn_rcv_tusd%cldes ), ' (', TRIM(sn_rcv_tusd%clcat ), ')' 
     372         WRITE(numout,*)' Transport associated to Stokes drift grid v = ', TRIM(sn_rcv_tvsd%cldes ), ' (', TRIM(sn_rcv_tvsd%clcat ), ')' 
     373         WRITE(numout,*)'      Bernouilli pressure head        = ', TRIM(sn_rcv_bhd%cldes   ), ' (', TRIM(sn_rcv_bhd%clcat  ), ')' 
     374         WRITE(numout,*)'Wave to ocean momentum flux and Net wave-supported stress = ', TRIM(sn_rcv_taw%cldes ), ' (', TRIM(sn_rcv_taw%clcat ), ')' 
     375         WRITE(numout,*)'      Surface current to waves        = ', TRIM(sn_snd_crtw%cldes  ), ' (', TRIM(sn_snd_crtw%clcat  ), ')' 
     376         WRITE(numout,*)'                      - referential   = ', sn_snd_crtw%clvref 
     377         WRITE(numout,*)'                      - orientation   = ', sn_snd_crtw%clvor 
     378         WRITE(numout,*)'                      - mesh          = ', sn_snd_crtw%clvgrd 
     379      ENDIF 
    333380      !                                   ! allocate sbccpl arrays 
    334381      IF( sbc_cpl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' ) 
     
    367414      IF(       TRIM( sn_rcv_tau%cldes ) == 'oce only' .OR. TRIM( sn_rcv_tau%cldes ) == 'oce and ice'  & 
    368415           .OR. TRIM( sn_rcv_tau%cldes ) == 'mixed oce-ice' ) THEN ! avoid working with the atmospheric fields if they are not coupled 
    369  
     416      ! 
    370417      IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' )   srcv(jpr_otx1:jpr_itz2)%nsgn = -1. 
    371418       
     
    608655         cpl_wper = .TRUE. 
    609656      ENDIF 
    610       srcv(jpr_wfreq)%clname = 'O_WFreq'     ! wave peak frequency  
    611       IF( TRIM(sn_rcv_wfreq%cldes ) == 'coupled' )  THEN 
    612          srcv(jpr_wfreq)%laction = .TRUE. 
    613          cpl_wfreq = .TRUE. 
    614       ENDIF 
    615657      srcv(jpr_wnum)%clname = 'O_WNum'       ! mean wave number 
    616658      IF( TRIM(sn_rcv_wnum%cldes ) == 'coupled' )  THEN 
     
    618660         cpl_wnum = .TRUE. 
    619661      ENDIF 
    620       srcv(jpr_tauwoc)%clname = 'O_TauOce'   ! stress fraction adsorbed by the wave 
    621       IF( TRIM(sn_rcv_tauwoc%cldes ) == 'coupled' )  THEN 
    622          srcv(jpr_tauwoc)%laction = .TRUE. 
    623          cpl_tauwoc = .TRUE. 
    624       ENDIF 
    625       srcv(jpr_tauwx)%clname = 'O_Tauwx'      ! ocean stress from wave in the x direction 
    626       srcv(jpr_tauwy)%clname = 'O_Tauwy'      ! ocean stress from wave in the y direction 
    627       IF( TRIM(sn_rcv_tauw%cldes ) == 'coupled' )  THEN 
    628          srcv(jpr_tauwx)%laction = .TRUE. 
    629          srcv(jpr_tauwy)%laction = .TRUE. 
    630          cpl_tauw = .TRUE. 
     662      srcv(jpr_wstrf)%clname = 'O_WStrf'     ! stress fraction adsorbed by the wave 
     663      IF( TRIM(sn_rcv_wstrf%cldes ) == 'coupled' )  THEN 
     664         srcv(jpr_wstrf)%laction = .TRUE. 
     665         cpl_wstrf = .TRUE. 
    631666      ENDIF 
    632667      srcv(jpr_wdrag)%clname = 'O_WDrag'     ! neutral surface drag coefficient 
     
    635670         cpl_wdrag = .TRUE. 
    636671      ENDIF 
    637       IF( srcv(jpr_tauwoc)%laction .AND. srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction ) & 
    638             CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', & 
    639                                      '(sn_rcv_tauwoc=coupled and sn_rcv_tauw=coupled)' ) 
     672      srcv(jpr_charn)%clname = 'O_Charn'     ! Chranock coefficient 
     673      IF( TRIM(sn_rcv_charn%cldes ) == 'coupled' )  THEN 
     674         srcv(jpr_charn)%laction = .TRUE. 
     675         cpl_charn = .TRUE. 
     676      ENDIF 
     677      srcv(jpr_bhd)%clname = 'O_Bhd'     ! Bernoulli head. waves' induced surface pressure 
     678      IF( TRIM(sn_rcv_bhd%cldes ) == 'coupled' )  THEN 
     679         srcv(jpr_bhd)%laction = .TRUE. 
     680         cpl_bhd = .TRUE. 
     681      ENDIF 
     682      srcv(jpr_tusd)%clname = 'O_Tusd'     ! zonal stokes transport 
     683      IF( TRIM(sn_rcv_tusd%cldes ) == 'coupled' )  THEN 
     684         srcv(jpr_tusd)%laction = .TRUE. 
     685         cpl_tusd = .TRUE. 
     686      ENDIF 
     687      srcv(jpr_tvsd)%clname = 'O_Tvsd'     ! meridional stokes tranmport 
     688      IF( TRIM(sn_rcv_tvsd%cldes ) == 'coupled' )  THEN 
     689         srcv(jpr_tvsd)%laction = .TRUE. 
     690         cpl_tvsd = .TRUE. 
     691      ENDIF 
     692 
     693      srcv(jpr_twox)%clname = 'O_Twox'     ! wave to ocean momentum flux in the u direction 
     694      srcv(jpr_twoy)%clname = 'O_Twoy'     ! wave to ocean momentum flux in the v direction 
     695      srcv(jpr_tawx)%clname = 'O_Tawx'     ! Net wave-supported stress in the u direction 
     696      srcv(jpr_tawy)%clname = 'O_Tawy'     ! Net wave-supported stress in the v direction 
     697      IF( TRIM(sn_rcv_taw%cldes ) == 'coupled' )  THEN 
     698         srcv(jpr_twox)%laction = .TRUE. 
     699         srcv(jpr_twoy)%laction = .TRUE. 
     700         srcv(jpr_tawx)%laction = .TRUE. 
     701         srcv(jpr_tawy)%laction = .TRUE. 
     702         cpl_taw = .TRUE. 
     703      ENDIF 
    640704      ! 
    641705      !                                                      ! ------------------------------- ! 
     
    698762         ! Change first letter to couple with atmosphere if already coupled OPA 
    699763         ! this is nedeed as each variable name used in the namcouple must be unique: 
    700          ! for example O_Runoff received by OPA from SAS and therefore O_Runoff received by SAS from the Atmosphere 
     764         ! for example O_Runoff received by OPA from SAS and therefore S_Runoff received by SAS from the Atmosphere 
    701765         DO jn = 1, jprcv 
    702766            IF( srcv(jn)%clname(1:1) == "O" ) srcv(jn)%clname = "S"//srcv(jn)%clname(2:LEN(srcv(jn)%clname)) 
     
    822886      END SELECT 
    823887 
     888      ! Initialise ice fractions from last coupling time to zero (needed by Met-Office) 
     889#if defined key_si3 || defined key_cice 
     890       a_i_last_couple(:,:,:) = 0._wp 
     891#endif 
    824892      !                                                      ! ------------------------- !  
    825893      !                                                      !      Ice Meltponds        !  
     
    10331101      !   initialisation of the coupler  ! 
    10341102      ! ================================ ! 
    1035  
    10361103      CALL cpl_define(jprcv, jpsnd, nn_cplmodel) 
    10371104       
     
    10461113      ENDIF 
    10471114      xcplmask(:,:,0) = 1. - SUM( xcplmask(:,:,1:nn_cplmodel), dim = 3 ) 
     1115      ! 
    10481116      ! 
    10491117   END SUBROUTINE sbc_cpl_init 
     
    11101178      REAL(wp) ::   zcdrag = 1.5e-3        ! drag coefficient 
    11111179      REAL(wp) ::   zzx, zzy               ! temporary variables 
    1112       REAL(wp), DIMENSION(jpi,jpj) ::   ztx, zty, zmsk, zemp, zqns, zqsr 
     1180      REAL(wp), DIMENSION(jpi,jpj) ::   ztx, zty, zmsk, zemp, zqns, zqsr, zcloud_fra 
    11131181      !!---------------------------------------------------------------------- 
    11141182      ! 
     
    11211189         IF( ncpl_qsr_freq /= 0) ncpl_qsr_freq = 86400 / ncpl_qsr_freq ! used by top 
    11221190          
     1191         IF ( ln_wave .AND. nn_components == 0 ) THEN 
     1192            ncpl_qsr_freq = 1; 
     1193            WRITE(numout,*) 'ncpl_qsr_freq is set to 1 when coupling NEMO with wave (without SAS) ' 
     1194         ENDIF 
    11231195      ENDIF 
    11241196      ! 
     
    11701242            !                               
    11711243            IF( srcv(jpr_otx1)%clgrid == 'T' ) THEN 
    1172                DO_2D( 0, 0, 0, 0 ) 
     1244               DO_2D( 0, 0, 0, 0 )                                        ! T ==> (U,V) 
    11731245                  frcv(jpr_otx1)%z3(ji,jj,1) = 0.5 * ( frcv(jpr_otx1)%z3(ji+1,jj  ,1) + frcv(jpr_otx1)%z3(ji,jj,1) ) 
    11741246                  frcv(jpr_oty1)%z3(ji,jj,1) = 0.5 * ( frcv(jpr_oty1)%z3(ji  ,jj+1,1) + frcv(jpr_oty1)%z3(ji,jj,1) ) 
     
    12241296         ENDIF 
    12251297      ENDIF 
    1226  
     1298!!$      !                                                      ! ========================= ! 
     1299!!$      SELECT CASE( TRIM( sn_rcv_clouds%cldes ) )             !       cloud fraction      ! 
     1300!!$      !                                                      ! ========================= ! 
     1301!!$      cloud_fra(:,:) = frcv(jpr_clfra)*z3(:,:,1) 
     1302!!$      END SELECT 
     1303!!$ 
     1304      zcloud_fra(:,:) = pp_cldf   ! should be real cloud fraction instead (as in the bulk) but needs to be read from atm. 
     1305      IF( ln_mixcpl ) THEN 
     1306         cloud_fra(:,:) = cloud_fra(:,:) * xcplmask(:,:,0) + zcloud_fra(:,:)* zmsk(:,:) 
     1307      ELSE 
     1308         cloud_fra(:,:) = zcloud_fra(:,:) 
     1309      ENDIF 
     1310      !                                                      ! ========================= ! 
    12271311      ! u(v)tau and taum will be modified by ice model 
    12281312      ! -> need to be reset before each call of the ice/fsbc       
     
    12831367         IF( srcv(jpr_hsig)%laction ) hsw(:,:) = frcv(jpr_hsig)%z3(:,:,1) 
    12841368      !  
    1285       !                                                      ! ========================= !   
    1286       !                                                      !    Wave peak frequency    !  
    1287       !                                                      ! ========================= !   
    1288          IF( srcv(jpr_wfreq)%laction ) wfreq(:,:) = frcv(jpr_wfreq)%z3(:,:,1) 
    1289       ! 
    12901369      !                                                      ! ========================= !  
    12911370      !                                                      !    Vertical mixing Qiao   ! 
     
    12941373 
    12951374         ! Calculate the 3D Stokes drift both in coupled and not fully uncoupled mode 
    1296          IF( srcv(jpr_sdrftx)%laction .OR. srcv(jpr_sdrfty)%laction .OR. srcv(jpr_wper)%laction & 
    1297                                       .OR. srcv(jpr_hsig)%laction   .OR. srcv(jpr_wfreq)%laction) THEN 
     1375         IF( srcv(jpr_sdrftx)%laction .OR. srcv(jpr_sdrfty)%laction .OR. & 
     1376             srcv(jpr_wper)%laction .OR. srcv(jpr_hsig)%laction )  THEN 
    12981377            CALL sbc_stokes( Kmm ) 
    12991378         ENDIF 
     
    13021381      !                                                      ! Stress adsorbed by waves  ! 
    13031382      !                                                      ! ========================= !  
    1304       IF( srcv(jpr_tauwoc)%laction .AND. ln_tauwoc ) tauoc_wave(:,:) = frcv(jpr_tauwoc)%z3(:,:,1) 
    1305  
    1306       !                                                      ! ========================= !   
    1307       !                                                      ! Stress component by waves !  
    1308       !                                                      ! ========================= !   
    1309       IF( srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction .AND. ln_tauw ) THEN 
    1310          tauw_x(:,:) = frcv(jpr_tauwx)%z3(:,:,1) 
    1311          tauw_y(:,:) = frcv(jpr_tauwy)%z3(:,:,1) 
    1312       ENDIF 
    1313  
     1383      IF( srcv(jpr_wstrf)%laction .AND. ln_tauoc )  tauoc_wave(:,:) = frcv(jpr_wstrf)%z3(:,:,1) 
     1384      ! 
    13141385      !                                                      ! ========================= !  
    13151386      !                                                      !   Wave drag coefficient   ! 
    13161387      !                                                      ! ========================= !  
    13171388      IF( srcv(jpr_wdrag)%laction .AND. ln_cdgw )   cdn_wave(:,:) = frcv(jpr_wdrag)%z3(:,:,1) 
    1318  
     1389      ! 
     1390      !                                                      ! ========================= ! 
     1391      !                                                      !   Chranock coefficient    ! 
     1392      !                                                      ! ========================= ! 
     1393      IF( srcv(jpr_charn)%laction .AND. ln_charn )  charn(:,:) = frcv(jpr_charn)%z3(:,:,1) 
     1394      ! 
     1395      !                                                      ! ========================= ! 
     1396      !                                                      ! net wave-supported stress ! 
     1397      !                                                      ! ========================= ! 
     1398      IF( srcv(jpr_tawx)%laction .AND. ln_taw )     tawx(:,:) = frcv(jpr_tawx)%z3(:,:,1) 
     1399      IF( srcv(jpr_tawy)%laction .AND. ln_taw )     tawy(:,:) = frcv(jpr_tawy)%z3(:,:,1) 
     1400      ! 
     1401      !                                                      ! ========================= ! 
     1402      !                                                      !wave to ocean momentum flux! 
     1403      !                                                      ! ========================= ! 
     1404      IF( srcv(jpr_twox)%laction .AND. ln_taw )     twox(:,:) = frcv(jpr_twox)%z3(:,:,1) 
     1405      IF( srcv(jpr_twoy)%laction .AND. ln_taw )     twoy(:,:) = frcv(jpr_twoy)%z3(:,:,1) 
     1406      !                                                       
     1407      !                                                      ! ========================= ! 
     1408      !                                                      !    wave TKE flux at sfc   ! 
     1409      !                                                      ! ========================= ! 
     1410      IF( srcv(jpr_phioc)%laction .AND. ln_phioc )     phioc(:,:) = frcv(jpr_phioc)%z3(:,:,1) 
     1411      ! 
     1412      !                                                      ! ========================= ! 
     1413      !                                                      !      Bernoulli head       ! 
     1414      !                                                      ! ========================= ! 
     1415      IF( srcv(jpr_bhd)%laction .AND. ln_bern_srfc )   bhd_wave(:,:) = frcv(jpr_bhd)%z3(:,:,1) 
     1416      ! 
     1417      !                                                      ! ========================= ! 
     1418      !                                                      !   Stokes transport u dir  ! 
     1419      !                                                      ! ========================= ! 
     1420      IF( srcv(jpr_tusd)%laction .AND. ln_breivikFV_2016 )    tusd(:,:) = frcv(jpr_tusd)%z3(:,:,1) 
     1421      ! 
     1422      !                                                      ! ========================= ! 
     1423      !                                                      !   Stokes transport v dir  ! 
     1424      !                                                      ! ========================= ! 
     1425      IF( srcv(jpr_tvsd)%laction .AND. ln_breivikFV_2016 )     tvsd(:,:) = frcv(jpr_tvsd)%z3(:,:,1) 
     1426      ! 
    13191427      !  Fields received by SAS when OASIS coupling 
    13201428      !  (arrays no more filled at sbcssm stage) 
     
    15491657            p_tauj(:,:) = frcv(jpr_ity1)%z3(:,:,1) 
    15501658         CASE( 'T' ) 
    1551             DO_2D( 0, 0, 0, 0 ) 
     1659            DO_2D( 0, 0, 0, 0 )                    ! T ==> (U,V) 
    15521660               ! take care of the land-sea mask to avoid "pollution" of coastal stress. p[uv]taui used in frazil and  rheology  
    15531661               zztmp1 = 0.5_wp * ( 2. - umask(ji,jj,1) ) * MAX( tmask(ji,jj,1),tmask(ji+1,jj  ,1) ) 
     
    16231731      ! 
    16241732      INTEGER  ::   ji, jj, jl   ! dummy loop index 
    1625       REAL(wp) ::   ztri         ! local scalar 
    16261733      REAL(wp), DIMENSION(jpi,jpj)     ::   zcptn, zcptrain, zcptsnw, ziceld, zmsk, zsnw 
    16271734      REAL(wp), DIMENSION(jpi,jpj)     ::   zemp_tot, zemp_ice, zemp_oce, ztprecip, zsprecip  , zevap_oce, zdevap_ice 
    16281735      REAL(wp), DIMENSION(jpi,jpj)     ::   zqns_tot, zqns_oce, zqsr_tot, zqsr_oce, zqprec_ice, zqemp_oce, zqemp_ice 
     1736      REAL(wp), DIMENSION(jpi,jpj)     ::   zevap_ice_total 
    16291737      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zqns_ice, zqsr_ice, zdqns_ice, zqevap_ice, zevap_ice, zqtr_ice_top, ztsu 
     1738      REAL(wp), DIMENSION(jpi,jpj)     ::   ztri 
    16301739      !!---------------------------------------------------------------------- 
    16311740      ! 
     
    16471756         ztprecip(:,:) =   frcv(jpr_rain)%z3(:,:,1) + zsprecip(:,:)  ! May need to ensure positive here 
    16481757         zemp_tot(:,:) =   frcv(jpr_tevp)%z3(:,:,1) - ztprecip(:,:) 
    1649          zemp_ice(:,:) = ( frcv(jpr_ievp)%z3(:,:,1) - frcv(jpr_snow)%z3(:,:,1) ) * picefr(:,:) 
    16501758      CASE( 'oce and ice'   )   ! received fields: jpr_sbpr, jpr_semp, jpr_oemp, jpr_ievp 
    16511759         zemp_tot(:,:) = ziceld(:,:) * frcv(jpr_oemp)%z3(:,:,1) + picefr(:,:) * frcv(jpr_sbpr)%z3(:,:,1) 
     
    16591767 
    16601768#if defined key_si3 
     1769 
     1770      ! --- evaporation over ice (kg/m2/s) --- ! 
     1771      IF (ln_scale_ice_flux) THEN ! typically met-office requirements 
     1772         IF (sn_rcv_emp%clcat == 'yes') THEN 
     1773            WHERE( a_i(:,:,:) > 1.e-10 )  ; zevap_ice(:,:,:) = frcv(jpr_ievp)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) 
     1774            ELSEWHERE                     ; zevap_ice(:,:,:) = 0._wp 
     1775            END WHERE 
     1776            WHERE( picefr(:,:) > 1.e-10 ) ; zevap_ice_total(:,:) = SUM( zevap_ice(:,:,:) * a_i(:,:,:), dim=3 ) / picefr(:,:) 
     1777            ELSEWHERE                     ; zevap_ice_total(:,:) = 0._wp 
     1778            END WHERE 
     1779         ELSE 
     1780            WHERE( picefr(:,:) > 1.e-10 ) ; zevap_ice(:,:,1) = frcv(jpr_ievp)%z3(:,:,1) * SUM( a_i_last_couple, dim=3 ) / picefr(:,:) 
     1781            ELSEWHERE                     ; zevap_ice(:,:,1) = 0._wp 
     1782            END WHERE 
     1783            zevap_ice_total(:,:) = zevap_ice(:,:,1) 
     1784            DO jl = 2, jpl 
     1785               zevap_ice(:,:,jl) = zevap_ice(:,:,1) 
     1786            ENDDO 
     1787         ENDIF 
     1788      ELSE 
     1789         IF (sn_rcv_emp%clcat == 'yes') THEN 
     1790            zevap_ice(:,:,1:jpl) = frcv(jpr_ievp)%z3(:,:,1:jpl) 
     1791            WHERE( picefr(:,:) > 1.e-10 ) ; zevap_ice_total(:,:) = SUM( zevap_ice(:,:,:) * a_i(:,:,:), dim=3 ) / picefr(:,:) 
     1792            ELSEWHERE                     ; zevap_ice_total(:,:) = 0._wp 
     1793            END WHERE 
     1794         ELSE 
     1795            zevap_ice(:,:,1) = frcv(jpr_ievp)%z3(:,:,1) 
     1796            zevap_ice_total(:,:) = zevap_ice(:,:,1) 
     1797            DO jl = 2, jpl 
     1798               zevap_ice(:,:,jl) = zevap_ice(:,:,1) 
     1799            ENDDO 
     1800         ENDIF 
     1801      ENDIF 
     1802 
     1803      IF ( TRIM( sn_rcv_emp%cldes ) == 'conservative' ) THEN 
     1804         ! For conservative case zemp_ice has not been defined yet. Do it now. 
     1805         zemp_ice(:,:) = zevap_ice_total(:,:) * picefr(:,:) - frcv(jpr_snow)%z3(:,:,1) * picefr(:,:) 
     1806      ENDIF 
     1807 
    16611808      ! zsnw = snow fraction over ice after wind blowing (=picefr if no blowing) 
    1662       zsnw(:,:) = 0._wp   ;   CALL ice_thd_snwblow( ziceld, zsnw ) 
     1809      zsnw(:,:) = 0._wp   ;   CALL ice_var_snwblow( ziceld, zsnw ) 
    16631810       
    16641811      ! --- evaporation minus precipitation corrected (because of wind blowing on snow) --- ! 
     
    16671814 
    16681815      ! --- evaporation over ocean (used later for qemp) --- ! 
    1669       zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) 
    1670  
    1671       ! --- evaporation over ice (kg/m2/s) --- ! 
    1672       DO jl=1,jpl 
    1673          IF(sn_rcv_emp%clcat == 'yes') THEN   ;   zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,jl) 
    1674          ELSE                                  ;   zevap_ice(:,:,jl) = frcv(jpr_ievp)%z3(:,:,1 )   ;   ENDIF 
    1675       ENDDO 
     1816      zevap_oce(:,:) = frcv(jpr_tevp)%z3(:,:,1) - zevap_ice_total(:,:) * picefr(:,:) 
    16761817 
    16771818      ! since the sensitivity of evap to temperature (devap/dT) is not prescribed by the atmosphere, we set it to 0 
     
    17511892!!      IF( srcv(jpr_rnf)%laction )   CALL iom_put( 'runoffs' , rnf(:,:) * tmask(:,:,1)                                 )  ! runoff 
    17521893!!      IF( srcv(jpr_isf)%laction )   CALL iom_put( 'iceshelf_cea', -fwfisf(:,:) * tmask(:,:,1)                         )  ! iceshelf 
    1753       IF( srcv(jpr_cal)%laction )   CALL iom_put( 'calving_cea' , frcv(jpr_cal)%z3(:,:,1) * tmask(:,:,1)                )  ! calving 
    1754       IF( srcv(jpr_icb)%laction )   CALL iom_put( 'iceberg_cea' , frcv(jpr_icb)%z3(:,:,1) * tmask(:,:,1)                )  ! icebergs 
    1755       IF( iom_use('snowpre') )      CALL iom_put( 'snowpre'     , sprecip(:,:)                                          )  ! Snow 
    1756       IF( iom_use('precip') )       CALL iom_put( 'precip'      , tprecip(:,:)                                          )  ! total  precipitation 
    1757       IF( iom_use('rain') )         CALL iom_put( 'rain'        , tprecip(:,:) - sprecip(:,:)                           )  ! liquid precipitation  
    1758       IF( iom_use('snow_ao_cea') )  CALL iom_put( 'snow_ao_cea' , sprecip(:,:) * ( 1._wp - zsnw(:,:) )                  )  ! Snow over ice-free ocean  (cell average) 
    1759       IF( iom_use('snow_ai_cea') )  CALL iom_put( 'snow_ai_cea' , sprecip(:,:) *           zsnw(:,:)                    )  ! Snow over sea-ice         (cell average) 
    1760       IF( iom_use('rain_ao_cea') )  CALL iom_put( 'rain_ao_cea' , ( tprecip(:,:) - sprecip(:,:) ) * picefr(:,:)         )  ! liquid precipitation over ocean (cell average) 
    1761       IF( iom_use('subl_ai_cea') )  CALL iom_put( 'subl_ai_cea' , frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * tmask(:,:,1) )  ! Sublimation over sea-ice (cell average) 
    1762       IF( iom_use('evap_ao_cea') )  CALL iom_put( 'evap_ao_cea' , ( frcv(jpr_tevp)%z3(:,:,1)  & 
    1763          &                                                        - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) ) * tmask(:,:,1) ) ! ice-free oce evap (cell average) 
     1894      IF( srcv(jpr_cal)%laction )    CALL iom_put( 'calving_cea' , frcv(jpr_cal)%z3(:,:,1) * tmask(:,:,1)                )  ! calving 
     1895      IF( srcv(jpr_icb)%laction )    CALL iom_put( 'iceberg_cea' , frcv(jpr_icb)%z3(:,:,1) * tmask(:,:,1)                )  ! icebergs 
     1896      IF( iom_use('snowpre') )       CALL iom_put( 'snowpre'     , sprecip(:,:)                                          )  ! Snow 
     1897      IF( iom_use('precip') )        CALL iom_put( 'precip'      , tprecip(:,:)                                          )  ! total  precipitation 
     1898      IF( iom_use('rain') )          CALL iom_put( 'rain'        , tprecip(:,:) - sprecip(:,:)                           )  ! liquid precipitation  
     1899      IF( iom_use('snow_ao_cea') )   CALL iom_put( 'snow_ao_cea' , sprecip(:,:) * ( 1._wp - zsnw(:,:) )                  )  ! Snow over ice-free ocean  (cell average) 
     1900      IF( iom_use('snow_ai_cea') )   CALL iom_put( 'snow_ai_cea' , sprecip(:,:) *           zsnw(:,:)                    )  ! Snow over sea-ice         (cell average) 
     1901      IF( iom_use('rain_ao_cea') )   CALL iom_put( 'rain_ao_cea' , ( tprecip(:,:) - sprecip(:,:) ) * picefr(:,:)         )  ! liquid precipitation over ocean (cell average) 
     1902      IF( iom_use('subl_ai_cea') )   CALL iom_put( 'subl_ai_cea' , frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) * tmask(:,:,1) )     ! Sublimation over sea-ice (cell average) 
     1903      IF( iom_use('evap_ao_cea') )   CALL iom_put( 'evap_ao_cea' , ( frcv(jpr_tevp)%z3(:,:,1)  & 
     1904         &                                                         - frcv(jpr_ievp)%z3(:,:,1) * picefr(:,:) ) * tmask(:,:,1) ) ! ice-free oce evap (cell average) 
    17641905      ! note: runoff output is done in sbcrnf (which includes icebergs too) and iceshelf output is done in sbcisf 
    17651906      ! 
     
    17691910      CASE( 'oce only' )         ! the required field is directly provided 
    17701911         zqns_tot(:,:) = frcv(jpr_qnsoce)%z3(:,:,1) 
     1912         ! For Met Office sea ice non-solar fluxes are already delt with by JULES so setting to zero 
     1913         ! here so the only flux is the ocean only one. 
     1914         zqns_ice(:,:,:) = 0._wp  
    17711915      CASE( 'conservative' )     ! the required fields are directly provided 
    17721916         zqns_tot(:,:) = frcv(jpr_qnsmix)%z3(:,:,1) 
     
    17981942               zqns_ice(:,:,jl) = frcv(jpr_qnsmix)%z3(:,:,jl)    & 
    17991943                  &             + frcv(jpr_dqnsdt)%z3(:,:,jl) * ( pist(:,:,jl) - ( ( rt0 + psst(:,:) ) * ziceld(:,:)   & 
    1800                   &                                                                + pist(:,:,jl) * picefr(:,:) ) ) 
     1944                  &                                             + pist(:,:,jl) * picefr(:,:) ) ) 
    18011945            END DO 
    18021946         ELSE 
     
    18041948               zqns_ice(:,:,jl) = frcv(jpr_qnsmix)%z3(:,:, 1)    & 
    18051949                  &             + frcv(jpr_dqnsdt)%z3(:,:, 1) * ( pist(:,:,jl) - ( ( rt0 + psst(:,:) ) * ziceld(:,:)   & 
    1806                   &                                                                + pist(:,:,jl) * picefr(:,:) ) ) 
     1950                  &                                             + pist(:,:,jl) * picefr(:,:) ) ) 
    18071951            END DO 
    18081952         ENDIF 
     
    19102054      CASE( 'oce only' ) 
    19112055         zqsr_tot(:,:  ) = MAX( 0._wp , frcv(jpr_qsroce)%z3(:,:,1) ) 
     2056         ! For Met Office sea ice solar fluxes are already delt with by JULES so setting to zero 
     2057         ! here so the only flux is the ocean only one. 
     2058         zqsr_ice(:,:,:) = 0._wp 
    19122059      CASE( 'conservative' ) 
    19132060         zqsr_tot(:,:  ) = frcv(jpr_qsrmix)%z3(:,:,1) 
     
    19952142            ENDDO 
    19962143         ENDIF 
     2144      CASE( 'none' )  
     2145         zdqns_ice(:,:,:) = 0._wp 
    19972146      END SELECT 
    19982147       
     
    20102159      !                                                      ! ========================= ! 
    20112160      CASE ('coupled') 
    2012          IF( ln_mixcpl ) THEN 
    2013             DO jl=1,jpl 
    2014                qml_ice(:,:,jl) = qml_ice(:,:,jl) * xcplmask(:,:,0) + frcv(jpr_topm)%z3(:,:,jl) * zmsk(:,:) 
    2015                qcn_ice(:,:,jl) = qcn_ice(:,:,jl) * xcplmask(:,:,0) + frcv(jpr_botm)%z3(:,:,jl) * zmsk(:,:) 
    2016             ENDDO 
     2161         IF (ln_scale_ice_flux) THEN 
     2162            WHERE( a_i(:,:,:) > 1.e-10_wp ) 
     2163               qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) 
     2164               qcn_ice(:,:,:) = frcv(jpr_botm)%z3(:,:,:) * a_i_last_couple(:,:,:) / a_i(:,:,:) 
     2165            ELSEWHERE 
     2166               qml_ice(:,:,:) = 0.0_wp 
     2167               qcn_ice(:,:,:) = 0.0_wp 
     2168            END WHERE 
    20172169         ELSE 
    20182170            qml_ice(:,:,:) = frcv(jpr_topm)%z3(:,:,:) 
     
    20252177      IF( .NOT.ln_cndflx ) THEN                              !==  No conduction flux as surface forcing  ==! 
    20262178         ! 
    2027          !                    ! ===> used prescribed cloud fraction representative for polar oceans in summer (0.81) 
    2028          ztri = 0.18 * ( 1.0 - cldf_ice ) + 0.35 * cldf_ice    ! surface transmission when hi>10cm (Grenfell Maykut 77) 
    2029          ! 
    2030          WHERE    ( phs(:,:,:) <= 0._wp .AND. phi(:,:,:) <  0.1_wp )       ! linear decrease from hi=0 to 10cm   
    2031             zqtr_ice_top(:,:,:) = qsr_ice(:,:,:) * ( ztri + ( 1._wp - ztri ) * ( 1._wp - phi(:,:,:) * 10._wp ) ) 
    2032          ELSEWHERE( phs(:,:,:) <= 0._wp .AND. phi(:,:,:) >= 0.1_wp )       ! constant (ztri) when hi>10cm 
    2033             zqtr_ice_top(:,:,:) = qsr_ice(:,:,:) * ztri 
    2034          ELSEWHERE                                                         ! zero when hs>0 
    2035             zqtr_ice_top(:,:,:) = 0._wp 
    2036          END WHERE 
     2179         IF( nn_qtrice == 0 ) THEN 
     2180            ! formulation derived from Grenfell and Maykut (1977), where transmission rate 
     2181            !    1) depends on cloudiness 
     2182            !       ! ===> used prescribed cloud fraction representative for polar oceans in summer (0.81) 
     2183            !       !      should be real cloud fraction instead (as in the bulk) but needs to be read from atm. 
     2184            !    2) is 0 when there is any snow 
     2185            !    3) tends to 1 for thin ice 
     2186            ztri(:,:) = 0.18 * ( 1.0 - cloud_fra(:,:) ) + 0.35 * cloud_fra(:,:)  ! surface transmission when hi>10cm 
     2187            DO jl = 1, jpl 
     2188               WHERE    ( phs(:,:,jl) <= 0._wp .AND. phi(:,:,jl) <  0.1_wp )       ! linear decrease from hi=0 to 10cm   
     2189                  zqtr_ice_top(:,:,jl) = zqsr_ice(:,:,jl) * ( ztri(:,:) + ( 1._wp - ztri(:,:) ) * ( 1._wp - phi(:,:,jl) * 10._wp ) ) 
     2190               ELSEWHERE( phs(:,:,jl) <= 0._wp .AND. phi(:,:,jl) >= 0.1_wp )       ! constant (ztri) when hi>10cm 
     2191                  zqtr_ice_top(:,:,jl) = zqsr_ice(:,:,jl) * ztri(:,:) 
     2192               ELSEWHERE                                                           ! zero when hs>0 
     2193                  zqtr_ice_top(:,:,jl) = 0._wp  
     2194               END WHERE 
     2195            ENDDO 
     2196         ELSEIF( nn_qtrice == 1 ) THEN 
     2197            ! formulation is derived from the thesis of M. Lebrun (2019). 
     2198            !    It represents the best fit using several sets of observations 
     2199            !    It comes with snow conductivities adapted to freezing/melting conditions (see icethd_zdf_bl99.F90) 
     2200            zqtr_ice_top(:,:,:) = 0.3_wp * zqsr_ice(:,:,:) 
     2201         ENDIF 
    20372202         !      
    20382203      ELSEIF( ln_cndflx .AND. .NOT.ln_cndemulate ) THEN      !==  conduction flux as surface forcing  ==! 
    20392204         ! 
    2040          !                    ! ===> here we must receive the qtr_ice_top array from the coupler 
    2041          !                           for now just assume zero (fully opaque ice) 
     2205         !          ! ===> here we must receive the qtr_ice_top array from the coupler 
     2206         !                 for now just assume zero (fully opaque ice) 
    20422207         zqtr_ice_top(:,:,:) = 0._wp 
    20432208         ! 
     
    20962261      ! 
    20972262      isec = ( kt - nit000 ) * NINT( rn_Dt )        ! date of exchanges 
     2263      info = OASIS_idle 
    20982264 
    20992265      zfr_l(:,:) = 1.- fr_i(:,:) 
     
    22342400      ENDIF 
    22352401 
     2402#if defined key_si3 || defined key_cice 
     2403      ! If this coupling was successful then save ice fraction for use between coupling points.  
     2404      ! This is needed for some calculations where the ice fraction at the last coupling point  
     2405      ! is needed.  
     2406      IF(  info == OASIS_Sent    .OR. info == OASIS_ToRest .OR. &  
     2407         & info == OASIS_SentOut .OR. info == OASIS_ToRestOut ) THEN  
     2408         IF ( sn_snd_thick%clcat == 'yes' ) THEN  
     2409           a_i_last_couple(:,:,1:jpl) = a_i(:,:,1:jpl) 
     2410         ENDIF 
     2411      ENDIF 
     2412#endif 
     2413 
    22362414      IF( ssnd(jps_fice1)%laction ) THEN 
    22372415         SELECT CASE( sn_snd_thick1%clcat ) 
     
    22972475            SELECT CASE( sn_snd_mpnd%clcat )   
    22982476            CASE( 'yes' )   
    2299                ztmp3(:,:,1:jpl) =  a_ip_frac(:,:,1:jpl) 
     2477               ztmp3(:,:,1:jpl) =  a_ip_eff(:,:,1:jpl) 
    23002478               ztmp4(:,:,1:jpl) =  h_ip(:,:,1:jpl)   
    23012479            CASE( 'no' )   
     
    23032481               ztmp4(:,:,:) = 0.0   
    23042482               DO jl=1,jpl   
    2305                  ztmp3(:,:,1) = ztmp3(:,:,1) + a_ip_frac(:,:,jpl)   
    2306                  ztmp4(:,:,1) = ztmp4(:,:,1) + h_ip(:,:,jpl)  
     2483                 ztmp3(:,:,1) = ztmp3(:,:,1) + a_ip_frac(:,:,jpl) 
     2484                 ztmp4(:,:,1) = ztmp4(:,:,1) + h_ip(:,:,jpl) 
    23072485               ENDDO   
    23082486            CASE default   ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_mpnd%clcat' )   
Note: See TracChangeset for help on using the changeset viewer.