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 11277 for branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90 – NEMO

Ignore:
Timestamp:
2019-07-17T15:29:15+02:00 (5 years ago)
Author:
kingr
Message:

Merged Juan's changes for running AMM15 woth wave coupling.
Corrected minor logic error to allow AMM7-uncoupled to reproduce earlier results.
Few line spacing changes to allow merging with OBS br and trunk rvn 5518.

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/UKMO/AMM15_v3_6_STABLE_package_collate/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r8058 r11277  
    2323   USE sbcapr 
    2424   USE sbcdcy          ! surface boundary condition: diurnal cycle 
     25   USE sbcwave         ! surface boundary condition: waves 
    2526   USE phycst          ! physical constants 
    2627#if defined key_lim3 
     
    4142   USE timing          ! Timing 
    4243   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
     44#if defined key_oasis3 || defined key_oasis3mct    
     45   USE mod_oasis                    ! OASIS3-MCT module    
     46#endif   
    4347   USE eosbn2 
    4448   USE sbcrnf   , ONLY : l_rnfcpl 
     
    105109   INTEGER, PARAMETER ::   jpr_e3t1st = 41            ! first T level thickness  
    106110   INTEGER, PARAMETER ::   jpr_fraqsr = 42            ! fraction of solar net radiation absorbed in the first ocean level 
    107    INTEGER, PARAMETER ::   jprcv      = 42            ! total number of fields received 
     111   INTEGER, PARAMETER ::   jpr_mslp   = 43            ! mean sea level pressure    
     112   INTEGER, PARAMETER ::   jpr_hsig   = 44            ! Hsig    
     113   INTEGER, PARAMETER ::   jpr_phioc  = 45            ! Wave=>ocean energy flux    
     114   INTEGER, PARAMETER ::   jpr_sdrftx = 46            ! Stokes drift on grid 1    
     115   INTEGER, PARAMETER ::   jpr_sdrfty = 47            ! Stokes drift on grid 2    
     116   INTEGER, PARAMETER ::   jpr_wper   = 48            ! Mean wave period   
     117   INTEGER, PARAMETER ::   jpr_wnum   = 49            ! Mean wavenumber   
     118   INTEGER, PARAMETER ::   jpr_tauoc  = 50            ! Stress fraction adsorbed by waves 
     119   INTEGER, PARAMETER ::   jpr_wdrag  = 51            ! Neutral surface drag coefficient   
     120   INTEGER, PARAMETER ::   jpr_wfreq  = 52            ! Wave peak frequency   
     121   INTEGER, PARAMETER ::   jpr_tauwx  = 53            ! x component of the ocean stress from waves  
     122   INTEGER, PARAMETER ::   jpr_tauwy  = 54            ! y component of the ocean stress from waves  
     123   INTEGER, PARAMETER ::   jprcv      = 54            ! total number of fields received 
    108124 
    109125   INTEGER, PARAMETER ::   jps_fice   =  1            ! ice fraction sent to the atmosphere 
     
    135151   INTEGER, PARAMETER ::   jps_e3t1st = 27            ! first level depth (vvl) 
    136152   INTEGER, PARAMETER ::   jps_fraqsr = 28            ! fraction of solar net radiation absorbed in the first ocean level 
    137    INTEGER, PARAMETER ::   jpsnd      = 28            ! total number of fields sended 
     153   INTEGER, PARAMETER ::   jps_ficet  = 29            ! total ice fraction     
     154   INTEGER, PARAMETER ::   jps_ocxw   = 30            ! currents on grid 1     
     155   INTEGER, PARAMETER ::   jps_ocyw   = 31            ! currents on grid 2   
     156   INTEGER, PARAMETER ::   jps_wlev   = 32            ! water level    
     157   INTEGER, PARAMETER ::   jpsnd      = 32            ! total number of fields sent 
    138158 
    139159   !                                                         !!** namelist namsbc_cpl ** 
     
    149169   ! Received from the atmosphere                     ! 
    150170   TYPE(FLD_C) ::   sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau, sn_rcv_dqnsdt, sn_rcv_qsr, sn_rcv_qns, sn_rcv_emp, sn_rcv_rnf 
    151    TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2                         
     171   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_mslp                              
     172   ! Send to waves    
     173   TYPE(FLD_C) ::   sn_snd_ifrac, sn_snd_crtw, sn_snd_wlev    
     174   ! Received from waves    
     175   TYPE(FLD_C) ::   sn_rcv_hsig,sn_rcv_phioc,sn_rcv_sdrft,sn_rcv_wper, &  
     176                    sn_rcv_wfreq,sn_rcv_wnum,sn_rcv_tauoc,sn_rcv_tauw, &  
     177                    sn_rcv_wdrag 
    152178   ! Other namelist parameters                        ! 
    153179   INTEGER     ::   nn_cplmodel            ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 
     
    161187 
    162188   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   albedo_oce_mix     ! ocean albedo sent to atmosphere (mix clear/overcast sky) 
    163  
     189     
    164190   INTEGER , ALLOCATABLE, SAVE, DIMENSION(    :) ::   nrcvinfo           ! OASIS info argument 
    165191 
     
    179205      !!             ***  FUNCTION sbc_cpl_alloc  *** 
    180206      !!---------------------------------------------------------------------- 
    181       INTEGER :: ierr(3) 
     207      INTEGER :: ierr(4) 
    182208      !!---------------------------------------------------------------------- 
    183209      ierr(:) = 0 
     
    188214      ALLOCATE( a_i(jpi,jpj,1) , STAT=ierr(2) )  ! used in sbcice_if.F90 (done here as there is no sbc_ice_if_init) 
    189215#endif 
    190       ALLOCATE( xcplmask(jpi,jpj,0:nn_cplmodel) , STAT=ierr(3) ) 
    191       ! 
     216!      ALLOCATE( xcplmask(jpi,jpj,0:nn_cplmodel) , STAT=ierr(3) ) 
     217      ! Hardwire three models as nn_cplmodel has not been read in from the namelist yet.    
     218      ALLOCATE( xcplmask(jpi,jpj,0:3) , STAT=ierr(3) ) 
     219      ! 
     220      IF( .NOT. ln_apr_dyn ) ALLOCATE( ssh_ib(jpi,jpj), ssh_ibb(jpi,jpj), apr(jpi, jpj), STAT=ierr(4) ) 
     221 
    192222      sbc_cpl_alloc = MAXVAL( ierr ) 
    193223      IF( lk_mpp            )   CALL mpp_sum ( sbc_cpl_alloc ) 
     
    216246      REAL(wp), POINTER, DIMENSION(:,:) ::   zacs, zaos 
    217247      !! 
    218       NAMELIST/namsbc_cpl/  sn_snd_temp, sn_snd_alb   , sn_snd_thick, sn_snd_crt   , sn_snd_co2,      & 
    219          &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau  , sn_rcv_dqnsdt, sn_rcv_qsr,      & 
    220          &                  sn_rcv_qns , sn_rcv_emp   , sn_rcv_rnf  , sn_rcv_cal   , sn_rcv_iceflx,   & 
    221          &                  sn_rcv_co2 , nn_cplmodel  , ln_usecplmask 
     248      NAMELIST/namsbc_cpl/  sn_snd_temp , sn_snd_alb  , sn_snd_thick , sn_snd_crt   , sn_snd_co2,      &    
     249         &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau   , sn_rcv_dqnsdt, sn_rcv_qsr,      &    
     250         &                  sn_snd_ifrac, sn_snd_crtw , sn_snd_wlev  , sn_rcv_hsig  , sn_rcv_phioc ,   &    
     251         &                  sn_rcv_sdrft, sn_rcv_wper  , sn_rcv_wnum  , sn_rcv_wfreq, sn_rcv_tauoc,    &  
     252         &                  sn_rcv_wdrag, sn_rcv_qns   , sn_rcv_emp   , sn_rcv_rnf  , sn_rcv_cal ,     &  
     253         &                  sn_rcv_iceflx, sn_rcv_co2   , sn_rcv_mslp , sn_rcv_tauw ,                  &  
     254         &                  nn_cplmodel, ln_usecplmask, ln_usernfmask 
    222255      !!--------------------------------------------------------------------- 
    223256      ! 
     
    260293         WRITE(numout,*)'      sea ice heat fluxes             = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')' 
    261294         WRITE(numout,*)'      atm co2                         = ', TRIM(sn_rcv_co2%cldes   ), ' (', TRIM(sn_rcv_co2%clcat   ), ')' 
     295         WRITE(numout,*)'      mean sea level pressure         = ', TRIM(sn_rcv_mslp%cldes  ), ' (', TRIM(sn_rcv_mslp%clcat  ), ')' 
     296         WRITE(numout,*)'      significant wave heigth         = ', TRIM(sn_rcv_hsig%cldes  ), ' (', TRIM(sn_rcv_hsig%clcat  ), ')'    
     297         WRITE(numout,*)'      wave to oce energy flux         = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')'    
     298         WRITE(numout,*)'      Surface Stokes drift u,v        = ', TRIM(sn_rcv_sdrft%cldes ), ' (', TRIM(sn_rcv_sdrft%clcat ), ')' 
     299         WRITE(numout,*)'      Mean wave period                = ', TRIM(sn_rcv_wper%cldes  ), ' (', TRIM(sn_rcv_wper%clcat  ), ')'    
     300         WRITE(numout,*)'      Mean wave number                = ', TRIM(sn_rcv_wnum%cldes  ), ' (', TRIM(sn_rcv_wnum%clcat  ), ')'    
     301         WRITE(numout,*)'      Wave peak frequency             = ', TRIM(sn_rcv_wfreq%cldes ), ' (', TRIM(sn_rcv_wfreq%clcat ), ')'    
     302         WRITE(numout,*)'      Stress frac adsorbed by waves   = ', TRIM(sn_rcv_tauoc%cldes ), ' (', TRIM(sn_rcv_tauoc%clcat ), ')'    
     303         WRITE(numout,*)'      Stress components by waves      = ', TRIM(sn_rcv_tauw%cldes  ), ' (', TRIM(sn_rcv_tauw%clcat  ), ')' 
     304         WRITE(numout,*)'      Neutral surf drag coefficient   = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')' 
    262305         WRITE(numout,*)'  sent fields (multiple ice categories)' 
    263306         WRITE(numout,*)'      surface temperature             = ', TRIM(sn_snd_temp%cldes  ), ' (', TRIM(sn_snd_temp%clcat  ), ')' 
    264307         WRITE(numout,*)'      albedo                          = ', TRIM(sn_snd_alb%cldes   ), ' (', TRIM(sn_snd_alb%clcat   ), ')' 
    265308         WRITE(numout,*)'      ice/snow thickness              = ', TRIM(sn_snd_thick%cldes ), ' (', TRIM(sn_snd_thick%clcat ), ')' 
     309         WRITE(numout,*)'      total ice fraction              = ', TRIM(sn_snd_ifrac%cldes ), ' (', TRIM(sn_snd_ifrac%clcat ), ')' 
    266310         WRITE(numout,*)'      surface current                 = ', TRIM(sn_snd_crt%cldes   ), ' (', TRIM(sn_snd_crt%clcat   ), ')' 
    267311         WRITE(numout,*)'                      - referential   = ', sn_snd_crt%clvref  
     
    269313         WRITE(numout,*)'                      - mesh          = ', sn_snd_crt%clvgrd 
    270314         WRITE(numout,*)'      oce co2 flux                    = ', TRIM(sn_snd_co2%cldes   ), ' (', TRIM(sn_snd_co2%clcat   ), ')' 
     315         WRITE(numout,*)'      water level                     = ', TRIM(sn_snd_wlev%cldes  ), ' (', TRIM(sn_snd_wlev%clcat  ), ')'    
     316         WRITE(numout,*)'      surface current to waves        = ', TRIM(sn_snd_crtw%cldes  ), ' (', TRIM(sn_snd_crtw%clcat  ), ')'    
     317         WRITE(numout,*)'                      - referential   = ', sn_snd_crtw%clvref    
     318         WRITE(numout,*)'                      - orientation   = ', sn_snd_crtw%clvor    
     319         WRITE(numout,*)'                      - mesh          = ', sn_snd_crtw%clvgrd 
    271320         WRITE(numout,*)'  nn_cplmodel                         = ', nn_cplmodel 
    272321         WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask 
     322         WRITE(numout,*)'  ln_usernfmask                       = ', ln_usernfmask 
    273323      ENDIF 
    274324 
    275325      !                                   ! allocate sbccpl arrays 
    276       IF( sbc_cpl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' ) 
     326      !IF( sbc_cpl_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'sbc_cpl_alloc : unable to allocate arrays' ) 
    277327      
    278328      ! ================================ ! 
     
    307357      !  
    308358      ! Vectors: change of sign at north fold ONLY if on the local grid 
     359      IF( TRIM( sn_rcv_tau%cldes ) == 'oce only' .OR. TRIM(sn_rcv_tau%cldes ) == 'oce and ice') THEN ! avoid working with the atmospheric fields if they are not coupled 
    309360      IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' )   srcv(jpr_otx1:jpr_itz2)%nsgn = -1. 
    310361       
     
    337388         srcv(jpr_otx2:jpr_otz2)%clgrid  = 'V'        !           and           V-point 
    338389         srcv(jpr_itx1:jpr_itz1)%clgrid  = 'F'        ! ice components given at F-point 
    339          srcv(jpr_otx1:jpr_otz2)%laction = .TRUE.     ! receive oce components on grid 1 & 2 
     390         !srcv(jpr_otx1:jpr_otz2)%laction = .TRUE.     ! receive oce components on grid 1 & 2   
     391         ! Currently needed for HadGEM3 - but shouldn't affect anyone else for the moment   
     392         srcv(jpr_otx1)%laction = .TRUE.    
     393         srcv(jpr_oty1)%laction = .TRUE.   
     394         !   
    340395         srcv(jpr_itx1:jpr_itz1)%laction = .TRUE.     ! receive ice components on grid 1 only 
    341396      CASE( 'T,I' )  
     
    374429         srcv(jpr_ity1)%clgrid = 'V'                  ! i.e. it is always at U- & V-points for i- & j-comp. resp. 
    375430      ENDIF 
     431      ENDIF 
    376432        
    377433      !                                                      ! ------------------------- ! 
     
    403459      IF( TRIM( sn_rcv_rnf%cldes ) == 'coupled' ) THEN 
    404460         srcv(jpr_rnf)%laction = .TRUE. 
    405          l_rnfcpl              = .TRUE.                      ! -> no need to read runoffs in sbcrnf 
     461         l_rnfcpl              = .NOT. ln_usernfmask         ! -> no need to read runoffs in sbcrnf 
    406462         ln_rnf                = nn_components /= jp_iam_sas ! -> force to go through sbcrnf if not sas 
    407463         IF(lwp) WRITE(numout,*) 
     
    470526      !                                                      ! ------------------------- ! 
    471527      srcv(jpr_co2 )%clname = 'O_AtmCO2'   ;   IF( TRIM(sn_rcv_co2%cldes   ) == 'coupled' )    srcv(jpr_co2 )%laction = .TRUE. 
     528 
     529      !                                                      ! ------------------------- !    
     530      !                                                      ! Mean Sea Level Pressure   !    
     531      !                                                      ! ------------------------- !    
     532      srcv(jpr_mslp)%clname = 'O_MSLP'  
     533      IF( TRIM(sn_rcv_mslp%cldes  ) == 'coupled' ) THEN  
     534         srcv(jpr_mslp)%laction = .TRUE.  
     535         cpl_mslp = .TRUE.  
     536      ENDIF 
     537 
    472538      !                                                      ! ------------------------- ! 
    473539      !                                                      !   topmelt and botmelt     !    
     
    483549         srcv(jpr_topm:jpr_botm)%laction = .TRUE. 
    484550      ENDIF 
     551      !                                                      ! ------------------------- !   
     552      !                                                      !      Wave breaking        !       
     553      !                                                      ! ------------------------- !    
     554      srcv(jpr_hsig)%clname  = 'O_Hsigwa'    ! significant wave height   
     555      IF( TRIM(sn_rcv_hsig%cldes  ) == 'coupled' )  THEN   
     556         srcv(jpr_hsig)%laction = .TRUE.   
     557         cpl_hsig = .TRUE.   
     558      ENDIF   
     559      srcv(jpr_phioc)%clname = 'O_PhiOce'    ! wave to ocean energy   
     560      IF( TRIM(sn_rcv_phioc%cldes ) == 'coupled' )  THEN   
     561         srcv(jpr_phioc)%laction = .TRUE.   
     562         cpl_phioc = .TRUE.   
     563      ENDIF   
     564      srcv(jpr_sdrftx)%clname = 'O_Sdrfx'    ! Stokes drift in the u direction   
     565      srcv(jpr_sdrfty)%clname = 'O_Sdrfy'    ! Stokes drift in the v direction   
     566      IF( TRIM(sn_rcv_sdrft%cldes ) == 'coupled' )  THEN 
     567         srcv(jpr_sdrftx)%laction = .TRUE.   
     568         srcv(jpr_sdrfty)%laction = .TRUE.   
     569         cpl_sdrft = .TRUE. 
     570      ENDIF   
     571      srcv(jpr_wper)%clname = 'O_WPer'       ! mean wave period   
     572      IF( TRIM(sn_rcv_wper%cldes  ) == 'coupled' )  THEN   
     573         srcv(jpr_wper)%laction = .TRUE.   
     574         cpl_wper = .TRUE.   
     575      ENDIF   
     576      srcv(jpr_wfreq)%clname = 'O_WFreq'     ! wave peak frequency   
     577      IF( TRIM(sn_rcv_wfreq%cldes ) == 'coupled' )  THEN   
     578         srcv(jpr_wfreq)%laction = .TRUE.   
     579         cpl_wfreq = .TRUE.   
     580      ENDIF 
     581      srcv(jpr_wnum)%clname = 'O_WNum'       ! mean wave number   
     582      IF( TRIM(sn_rcv_wnum%cldes ) == 'coupled' )  THEN   
     583         srcv(jpr_wnum)%laction = .TRUE.   
     584         cpl_wnum = .TRUE.   
     585      ENDIF   
     586      srcv(jpr_tauoc)%clname = 'O_TauOce'     ! stress fraction adsorbed by the wave   
     587      IF( TRIM(sn_rcv_tauoc%cldes ) == 'coupled' )  THEN   
     588         srcv(jpr_tauoc)%laction = .TRUE.   
     589         cpl_tauoc = .TRUE.   
     590      ENDIF   
     591      srcv(jpr_tauwx)%clname = 'O_Tauwx'      ! ocean stress from wave in the x direction  
     592      srcv(jpr_tauwy)%clname = 'O_Tauwy'      ! ocean stress from wave in the y direction  
     593      IF( TRIM(sn_rcv_tauw%cldes ) == 'coupled' )  THEN   
     594         srcv(jpr_tauwx)%laction = .TRUE.   
     595         srcv(jpr_tauwy)%laction = .TRUE.   
     596         cpl_tauw = .TRUE. 
     597      ENDIF   
     598      srcv(jpr_wdrag)%clname = 'O_WDrag'     ! neutral surface drag coefficient   
     599      IF( TRIM(sn_rcv_wdrag%cldes ) == 'coupled' )  THEN   
     600         srcv(jpr_wdrag)%laction = .TRUE.   
     601         cpl_wdrag = .TRUE.   
     602      ENDIF   
     603      !  
     604      IF( srcv(jpr_tauoc)%laction .AND. srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction ) &  
     605            CALL ctl_stop( 'More than one method for modifying the ocean stress has been selected ', &  
     606                                     '(sn_rcv_tauoc=coupled and sn_rcv_tauw=coupled)' )  
     607      ! 
     608      !    
    485609      !                                                      ! ------------------------------- ! 
    486610      !                                                      !   OPA-SAS coupling - rcv by opa !    
     
    637761      !                                                      ! ------------------------- ! 
    638762      ssnd(jps_fice)%clname = 'OIceFrc' 
     763      ssnd(jps_ficet)%clname = 'OIceFrcT' 
    639764      ssnd(jps_hice)%clname = 'OIceTck' 
    640765      ssnd(jps_hsnw)%clname = 'OSnwTck' 
     
    645770      ENDIF 
    646771       
     772      IF (TRIM( sn_snd_ifrac%cldes )  == 'coupled') ssnd(jps_ficet)%laction = .TRUE. 
     773 
    647774      SELECT CASE ( TRIM( sn_snd_thick%cldes ) ) 
    648775      CASE( 'none'         )       ! nothing to do 
     
    665792      ssnd(jps_ocy1)%clname = 'O_OCury1'   ;   ssnd(jps_ivy1)%clname = 'O_IVely1' 
    666793      ssnd(jps_ocz1)%clname = 'O_OCurz1'   ;   ssnd(jps_ivz1)%clname = 'O_IVelz1' 
     794      ssnd(jps_ocxw)%clname = 'O_OCurxw'    
     795      ssnd(jps_ocyw)%clname = 'O_OCuryw' 
    667796      ! 
    668797      ssnd(jps_ocx1:jps_ivz1)%nsgn = -1.   ! vectors: change of the sign at the north fold 
     
    685814      END SELECT 
    686815 
     816      ssnd(jps_ocxw:jps_ocyw)%nsgn = -1.   ! vectors: change of the sign at the north fold    
     817               
     818      IF( sn_snd_crtw%clvgrd == 'U,V' ) THEN    
     819         ssnd(jps_ocxw)%clgrid = 'U' ; ssnd(jps_ocyw)%clgrid = 'V'    
     820      ELSE IF( sn_snd_crtw%clvgrd /= 'T' ) THEN    
     821         CALL ctl_stop( 'sn_snd_crtw%clvgrd must be equal to T' )    
     822      ENDIF    
     823      IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) ssnd(jps_ocxw:jps_ocyw)%nsgn = 1.    
     824      SELECT CASE( TRIM( sn_snd_crtw%cldes ) )    
     825         CASE( 'none'                 )   ; ssnd(jps_ocxw:jps_ocyw)%laction = .FALSE.    
     826         CASE( 'oce only'             )   ; ssnd(jps_ocxw:jps_ocyw)%laction = .TRUE.    
     827         CASE( 'weighted oce and ice' )   !   nothing to do    
     828         CASE( 'mixed oce-ice'        )   ; ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE.    
     829         CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_crtw%cldes' )    
     830      END SELECT 
     831 
    687832      !                                                      ! ------------------------- ! 
    688833      !                                                      !          CO2 flux         ! 
    689834      !                                                      ! ------------------------- ! 
    690835      ssnd(jps_co2)%clname = 'O_CO2FLX' ;  IF( TRIM(sn_snd_co2%cldes) == 'coupled' )    ssnd(jps_co2 )%laction = .TRUE. 
     836 
     837      !                                                      ! ------------------------- !    
     838      !                                                      !     Sea surface height    !    
     839      !                                                      ! ------------------------- !    
     840      ssnd(jps_wlev)%clname = 'O_Wlevel' ;  IF( TRIM(sn_snd_wlev%cldes) == 'coupled' )   ssnd(jps_wlev)%laction = .TRUE. 
    691841 
    692842      !                                                      ! ------------------------------- ! 
     
    783933      IF( ln_dm2dc .AND. ln_cpl .AND. ncpl_qsr_freq /= 86400 )   & 
    784934         &   CALL ctl_stop( 'sbc_cpl_init: diurnal cycle reconstruction (ln_dm2dc) needs daily couping for solar radiation' ) 
    785       ncpl_qsr_freq = 86400 / ncpl_qsr_freq 
     935      IF( ln_dm2dc .AND. ln_cpl ) ncpl_qsr_freq = 86400 / ncpl_qsr_freq 
    786936 
    787937      CALL wrk_dealloc( jpi,jpj, zacs, zaos ) 
     
    837987      !!                        emp          upward mass flux [evap. - precip. (- runoffs) (- calving)] (ocean only case) 
    838988      !!---------------------------------------------------------------------- 
     989      USE sbcflx ,  ONLY : ln_shelf_flx  
     990      USE sbcssm ,  ONLY : sbc_ssm_cpl  
     991      USE lib_fortran     ! distributed memory computing library 
     992 
    839993      INTEGER, INTENT(in)           ::   kt          ! ocean model time step index 
    840994      INTEGER, INTENT(in)           ::   k_fsbc      ! frequency of sbc (-> ice model) computation  
     
    845999      INTEGER  ::   ji, jj, jn             ! dummy loop indices 
    8461000      INTEGER  ::   isec                   ! number of seconds since nit000 (assuming rdttra did not change since nit000) 
     1001      INTEGER  ::   ikchoix 
    8471002      REAL(wp) ::   zcumulneg, zcumulpos   ! temporary scalars      
    8481003      REAL(wp) ::   zcoef                  ! temporary scalar 
     
    8501005      REAL(wp) ::   zcdrag = 1.5e-3        ! drag coefficient 
    8511006      REAL(wp) ::   zzx, zzy               ! temporary variables 
    852       REAL(wp), POINTER, DIMENSION(:,:) ::   ztx, zty, zmsk, zemp, zqns, zqsr 
     1007      REAL(wp), POINTER, DIMENSION(:,:) ::   ztx, zty, ztx2, zty2, zmsk, zemp, zqns, zqsr 
    8531008      !!---------------------------------------------------------------------- 
    8541009      ! 
    8551010      IF( nn_timing == 1 )  CALL timing_start('sbc_cpl_rcv') 
    8561011      ! 
    857       CALL wrk_alloc( jpi,jpj, ztx, zty, zmsk, zemp, zqns, zqsr ) 
     1012      CALL wrk_alloc( jpi,jpj, ztx, zty, ztx2, zty2, zmsk, zemp, zqns, zqsr ) 
    8581013      ! 
    8591014      IF( ln_mixcpl )   zmsk(:,:) = 1. - xcplmask(:,:,0) 
     
    8931048            IF( TRIM( sn_rcv_tau%clvor ) == 'eastward-northward' ) THEN   ! 2 components oriented along the local grid 
    8941049               !                                                       ! (geographical to local grid -> rotate the components) 
    895                CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->i', ztx )    
    896                IF( srcv(jpr_otx2)%laction ) THEN 
    897                   CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty )    
    898                ELSE   
    899                   CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty )   
     1050               IF( srcv(jpr_otx1)%clgrid == 'U' .AND. (.NOT. srcv(jpr_otx2)%laction) ) THEN    
     1051                  ! Temporary code for HadGEM3 - will be removed eventually.    
     1052                  ! Only applies when we have only taux on U grid and tauy on V grid    
     1053                  DO jj=2,jpjm1    
     1054                     DO ji=2,jpim1    
     1055                        ztx(ji,jj)=0.25*vmask(ji,jj,1) &    
     1056                           *(frcv(jpr_otx1)%z3(ji,jj,1)+frcv(jpr_otx1)%z3(ji-1,jj,1)    &    
     1057                           +frcv(jpr_otx1)%z3(ji,jj+1,1)+frcv(jpr_otx1)%z3(ji-1,jj+1,1))    
     1058                        zty(ji,jj)=0.25*umask(ji,jj,1) &    
     1059                           *(frcv(jpr_oty1)%z3(ji,jj,1)+frcv(jpr_oty1)%z3(ji+1,jj,1)    &    
     1060                           +frcv(jpr_oty1)%z3(ji,jj-1,1)+frcv(jpr_oty1)%z3(ji+1,jj-1,1))    
     1061                     ENDDO    
     1062                  ENDDO    
     1063                               
     1064                  ikchoix = 1    
     1065                  CALL repcmo(frcv(jpr_otx1)%z3(:,:,1),zty,ztx,frcv(jpr_oty1)%z3(:,:,1),ztx2,zty2,ikchoix)    
     1066                  CALL lbc_lnk (ztx2,'U', -1. )    
     1067                  CALL lbc_lnk (zty2,'V', -1. )    
     1068                  frcv(jpr_otx1)%z3(:,:,1)=ztx2(:,:)    
     1069                  frcv(jpr_oty1)%z3(:,:,1)=zty2(:,:)    
     1070               ELSE    
     1071                  CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->i', ztx )       
     1072                  frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:)      ! overwrite 1st component on the 1st grid    
     1073                  IF( srcv(jpr_otx2)%laction ) THEN    
     1074                     CALL rot_rep( frcv(jpr_otx2)%z3(:,:,1), frcv(jpr_oty2)%z3(:,:,1), srcv(jpr_otx2)%clgrid, 'en->j', zty )       
     1075                  ELSE    
     1076                     CALL rot_rep( frcv(jpr_otx1)%z3(:,:,1), frcv(jpr_oty1)%z3(:,:,1), srcv(jpr_otx1)%clgrid, 'en->j', zty )     
     1077                  ENDIF    
     1078                  frcv(jpr_oty1)%z3(:,:,1) = zty(:,:)      ! overwrite 2nd component on the 2nd grid 
    9001079               ENDIF 
    901                frcv(jpr_otx1)%z3(:,:,1) = ztx(:,:)      ! overwrite 1st component on the 1st grid 
    902                frcv(jpr_oty1)%z3(:,:,1) = zty(:,:)      ! overwrite 2nd component on the 2nd grid 
    9031080            ENDIF 
    9041081            !                               
     
    9191096      ELSE                                                   !   No dynamical coupling   ! 
    9201097         !                                                   ! ========================= ! 
     1098         ! it is possible that the momentum is calculated from the winds (ln_shelf_flx) and a coupled drag coefficient  
     1099         IF( srcv(jpr_wdrag)%laction .AND. ln_shelf_flx .AND. ln_cdgw .AND. nn_drag == jp_std ) THEN  
     1100            DO jj = 1, jpjm1  
     1101               DO ji = 1, jpim1  
     1102                  ! here utau and vtau should contain the wind components as read from the forcing files,  
     1103                  ! and the wind module should already be properly calculated  
     1104                  frcv(jpr_otx1)%z3(ji,jj,1) = zrhoa * 0.5 * ( frcv(jpr_wdrag)%z3(ji,jj,1) + frcv(jpr_wdrag)%z3(ji+1,jj,1) ) * &  
     1105                                                                         utau(ji,jj) * 0.5 * ( wndm(ji,jj) + wndm(ji+1,jj) )  
     1106                  frcv(jpr_oty1)%z3(ji,jj,1) = zrhoa * 0.5 * ( frcv(jpr_wdrag)%z3(ji,jj,1) + frcv(jpr_wdrag)%z3(ji,jj+1,1) ) * &  
     1107                                                                         vtau(ji,jj) * 0.5 * ( wndm(ji,jj) + wndm(ji,jj+1) )  
     1108                  utau(ji,jj) = frcv(jpr_otx1)%z3(ji,jj,1)  
     1109                  vtau(ji,jj) = frcv(jpr_oty1)%z3(ji,jj,1)  
     1110               END DO  
     1111            END DO  
     1112            CALL lbc_lnk_multi( frcv(jpr_otx1)%z3(:,:,1), 'U', -1. , frcv(jpr_oty1)%z3(:,:,1), 'V', -1. , &  
     1113                                                             utau(:,:), 'U', -1. , vtau(:,:), 'V',  -1. )  
     1114            llnewtx = .TRUE.  
     1115         ELSE 
    9211116         frcv(jpr_otx1)%z3(:,:,1) = 0.e0                               ! here simply set to zero  
    9221117         frcv(jpr_oty1)%z3(:,:,1) = 0.e0                               ! an external read in a file can be added instead 
    9231118         llnewtx = .TRUE. 
     1119         ENDIF 
    9241120         ! 
    9251121      ENDIF 
     
    9411137            END DO 
    9421138            CALL lbc_lnk( frcv(jpr_taum)%z3(:,:,1), 'T', 1. ) 
     1139            IF( .NOT. srcv(jpr_otx1)%laction .AND. srcv(jpr_wdrag)%laction .AND. &  
     1140                                ln_shelf_flx .AND. ln_cdgw .AND. nn_drag == jp_std ) &  
     1141               taum(:,:) = frcv(jpr_taum)%z3(:,:,1) 
    9431142            llnewtau = .TRUE. 
    9441143         ELSE 
     
    9551154      !                                                      ! ========================= ! 
    9561155      !                                                      !      10 m wind speed      !   (wndm) 
     1156      !                                                      !   include wave drag coef  !   (wndm) 
    9571157      !                                                      ! ========================= ! 
    9581158      ! 
     
    9651165!CDIR NOVERRCHK 
    9661166               DO ji = 1, jpi  
     1167                  IF( ln_shelf_flx ) THEN   ! the 10 wind module is properly calculated before if ln_shelf_flx  
     1168                     frcv(jpr_w10m)%z3(ji,jj,1) = wndm(ji,jj)  
     1169                  ELSE 
    9671170                  frcv(jpr_w10m)%z3(ji,jj,1) = SQRT( frcv(jpr_taum)%z3(ji,jj,1) * zcoef ) 
     1171                  ENDIF 
    9681172               END DO 
    9691173            END DO 
     
    9751179      IF( MOD( kt-1, k_fsbc ) == 0 ) THEN 
    9761180         ! 
     1181         ! if ln_wavcpl, the fields already contain the right information from forcing even if not ln_mixcpl 
    9771182         IF( ln_mixcpl ) THEN 
    978             utau(:,:) = utau(:,:) * xcplmask(:,:,0) + frcv(jpr_otx1)%z3(:,:,1) * zmsk(:,:) 
    979             vtau(:,:) = vtau(:,:) * xcplmask(:,:,0) + frcv(jpr_oty1)%z3(:,:,1) * zmsk(:,:) 
    980             taum(:,:) = taum(:,:) * xcplmask(:,:,0) + frcv(jpr_taum)%z3(:,:,1) * zmsk(:,:) 
    981             wndm(:,:) = wndm(:,:) * xcplmask(:,:,0) + frcv(jpr_w10m)%z3(:,:,1) * zmsk(:,:) 
    982          ELSE 
     1183            IF( srcv(jpr_otx1)%laction ) THEN  
     1184               utau(:,:) = utau(:,:) * xcplmask(:,:,0) + frcv(jpr_otx1)%z3(:,:,1) * zmsk(:,:)  
     1185               vtau(:,:) = vtau(:,:) * xcplmask(:,:,0) + frcv(jpr_oty1)%z3(:,:,1) * zmsk(:,:)  
     1186            ENDIF  
     1187            IF( srcv(jpr_taum)%laction )   &  
     1188               taum(:,:) = taum(:,:) * xcplmask(:,:,0) + frcv(jpr_taum)%z3(:,:,1) * zmsk(:,:)  
     1189            IF( srcv(jpr_w10m)%laction )   &  
     1190               wndm(:,:) = wndm(:,:) * xcplmask(:,:,0) + frcv(jpr_w10m)%z3(:,:,1) * zmsk(:,:)  
     1191         ELSE IF( ll_purecpl ) THEN 
    9831192            utau(:,:) = frcv(jpr_otx1)%z3(:,:,1) 
    9841193            vtau(:,:) = frcv(jpr_oty1)%z3(:,:,1) 
     
    9961205      IF( srcv(jpr_co2)%laction )   atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1) 
    9971206#endif 
     1207      !    
     1208      !                                                      ! ========================= !    
     1209      !                                                      ! Mean Sea Level Pressure   !   (taum)    
     1210      !                                                      ! ========================= !    
     1211      !    
     1212      IF( srcv(jpr_mslp)%laction ) THEN                    ! UKMO SHELF effect of atmospheric pressure on SSH    
     1213          IF( kt /= nit000 )   ssh_ibb(:,:) = ssh_ib(:,:)    !* Swap of ssh_ib fields    
     1214        
     1215          !                                                  !* update the reference atmospheric pressure (if necessary)  
     1216          IF( ln_ref_apr )  rn_pref = glob_sum( frcv(jpr_mslp)%z3(:,:,1) * e1e2t(:,:) ) / tarea  
     1217         
     1218          ssh_ib(:,:) = - ( frcv(jpr_mslp)%z3(:,:,1) - rn_pref ) * r1_grau    ! equivalent ssh (inverse barometer)    
     1219          apr   (:,:) =     frcv(jpr_mslp)%z3(:,:,1) !atmospheric pressure    
     1220          !  
     1221          CALL iom_put( "ssh_ib", ssh_ib )                                    !* output the inverse barometer ssh 
     1222        
     1223          !                                         ! ---------------------------------------- !  
     1224          IF( kt == nit000 ) THEN                   !   set the forcing field at nit000 - 1    !  
     1225             !                                      ! ---------------------------------------- !  
     1226             !* Restart: read in restart file  
     1227             IF( ln_rstart .AND. iom_varid( numror, 'ssh_ibb', ldstop = .FALSE. ) > 0 ) THEN  
     1228                IF(lwp) WRITE(numout,*) 'sbc_cpl:   ssh_ibb read in the restart file'  
     1229                CALL iom_get( numror, jpdom_autoglo, 'ssh_ibb', ssh_ibb )   ! before inv. barometer ssh  
     1230             ELSE                                         !* no restart: set from nit000 values  
     1231                IF(lwp) WRITE(numout,*) 'sbc_cpl:   ssh_ibb set to nit000 values'  
     1232                ssh_ibb(:,:) = ssh_ib(:,:)  
     1233             ENDIF  
     1234          ENDIF  
     1235          !                                         ! ---------------------------------------- !  
     1236          IF( lrst_oce ) THEN                       !      Write in the ocean restart file     !  
     1237             !                                      ! ---------------------------------------- !  
     1238             IF(lwp) WRITE(numout,*)  
     1239             IF(lwp) WRITE(numout,*) 'sbc_cpl : ssh_ib written in ocean restart file at it= ', kt,' date= ', ndastp  
     1240             IF(lwp) WRITE(numout,*) '~~~~'  
     1241             CALL iom_rstput( kt, nitrst, numrow, 'ssh_ibb' , ssh_ib )  
     1242          ENDIF  
     1243         
     1244          ! Update mean ssh  
     1245          CALL sbc_ssm_cpl( kt ) 
     1246      END IF    
     1247      !   
     1248      IF( ln_sdw ) THEN  ! Stokes Drift correction activated   
     1249      !                                                      ! ========================= !    
     1250      !                                                      !       Stokes drift u,v    !   
     1251      !                                                      ! ========================= !    
     1252      IF( srcv(jpr_sdrftx)%laction .AND. srcv(jpr_sdrfty)%laction ) THEN  
     1253                                     ut0sd(:,:) = frcv(jpr_sdrftx)%z3(:,:,1)   
     1254                                     vt0sd(:,:) = frcv(jpr_sdrfty)%z3(:,:,1)   
     1255      ENDIF 
     1256      !   
     1257      !                                                      ! ========================= !    
     1258      !                                                      !      Wave mean period     !   
     1259      !                                                      ! ========================= !    
     1260         IF( srcv(jpr_wper)%laction ) wmp(:,:) = frcv(jpr_wper)%z3(:,:,1)   
     1261      !   
     1262      !                                                      ! ========================= !    
     1263      !                                                      !  Significant wave height  !   
     1264      !                                                      ! ========================= !    
     1265         IF( srcv(jpr_hsig)%laction ) hsw(:,:) = frcv(jpr_hsig)%z3(:,:,1)   
     1266      !   
     1267      !                                                      ! ========================= !    
     1268      !                                                      !    Wave peak frequency    !   
     1269      !                                                      ! ========================= !    
     1270         IF( srcv(jpr_wfreq)%laction ) wfreq(:,:) = frcv(jpr_wfreq)%z3(:,:,1)   
     1271      !  
     1272      !                                                      ! ========================= ! 
     1273      !                                                      !    Vertical mixing Qiao   !   
     1274      !                                                      ! ========================= !    
     1275         IF( srcv(jpr_wnum)%laction .AND. ln_zdfqiao ) wnum(:,:) = frcv(jpr_wnum)%z3(:,:,1)   
     1276        
     1277         ! Calculate the 3D Stokes drift both in coupled and not fully uncoupled mode   
     1278         IF( (srcv(jpr_sdrftx)%laction .AND. srcv(jpr_sdrfty)%laction) .OR. srcv(jpr_wper)%laction &   
     1279                                        .OR. srcv(jpr_hsig)%laction   .OR. srcv(jpr_wfreq)%laction) & 
     1280            CALL sbc_stokes()   
     1281      ENDIF   
     1282      !                                                      ! ========================= !    
     1283      !                                                      ! Stress adsorbed by waves  !   
     1284      !                                                      ! ========================= !    
     1285      IF( srcv(jpr_tauoc)%laction .AND. ln_tauoc ) THEN  
     1286         tauoc_wave(:,:) = frcv(jpr_tauoc)%z3(:,:,1)  
     1287         ! cap the value of tauoc  
     1288         WHERE(tauoc_wave <   0.0 ) tauoc_wave = 1.0  
     1289         WHERE(tauoc_wave > 100.0 ) tauoc_wave = 1.0  
     1290      ENDIF  
     1291      !                                                      ! ========================= !    
     1292      !                                                      ! Stress component by waves !   
     1293      !                                                      ! ========================= !    
     1294      IF( srcv(jpr_tauwx)%laction .AND. srcv(jpr_tauwy)%laction .AND. ln_tauw ) THEN  
     1295         tauw_x(:,:) = frcv(jpr_tauwx)%z3(:,:,1)  
     1296         tauw_y(:,:) = frcv(jpr_tauwy)%z3(:,:,1)  
     1297         ! cap the value of tauoc  
     1298         WHERE(tauw_x < -100.0 ) tauw_x = 0.0  
     1299         WHERE(tauw_x >  100.0 ) tauw_x = 0.0  
     1300         WHERE(tauw_y < -100.0 ) tauw_y = 0.0  
     1301         WHERE(tauw_y >  100.0 ) tauw_y = 0.0  
     1302      ENDIF 
     1303        
     1304      !                                                      ! ========================= !    
     1305      !                                                      !   Wave to ocean energy    ! 
     1306      !                                                      ! ========================= !    
     1307      IF( srcv(jpr_phioc)%laction .AND. ln_phioc ) THEN  
     1308         rn_crban(:,:) = 29.0 * frcv(jpr_phioc)%z3(:,:,1)  
     1309         WHERE( rn_crban <    0.0 ) rn_crban = 0.0  
     1310         WHERE( rn_crban > 1000.0 ) rn_crban = 1000.0  
     1311      ENDIF  
    9981312 
    9991313      !  Fields received by SAS when OASIS coupling 
     
    10671381               CALL ctl_stop( 'sbc_cpl_rcv: wrong definition of sn_rcv_emp%cldes' ) 
    10681382            END SELECT 
    1069          ELSE 
     1383         ELSE IF( ll_purecpl ) THEN 
    10701384            zemp(:,:) = 0._wp 
    10711385         ENDIF 
     
    10751389         IF( srcv(jpr_cal)%laction )     zemp(:,:) = zemp(:,:) - frcv(jpr_cal)%z3(:,:,1) 
    10761390          
    1077          IF( ln_mixcpl ) THEN   ;   emp(:,:) = emp(:,:) * xcplmask(:,:,0) + zemp(:,:) * zmsk(:,:) 
    1078          ELSE                   ;   emp(:,:) =                              zemp(:,:) 
     1391         IF( ln_mixcpl .AND. ( srcv(jpr_oemp)%laction .OR. srcv(jpr_rain)%laction )) THEN  
     1392                                         emp(:,:) = emp(:,:) * xcplmask(:,:,0) + zemp(:,:) * zmsk(:,:)  
     1393         ELSE IF( ll_purecpl ) THEN  ;   emp(:,:) = zemp(:,:) 
    10791394         ENDIF 
    10801395         ! 
     
    10911406            ENDIF 
    10921407         ENDIF 
    1093          IF( ln_mixcpl ) THEN   ;   qns(:,:) = qns(:,:) * xcplmask(:,:,0) + zqns(:,:) * zmsk(:,:) 
    1094          ELSE                   ;   qns(:,:) =                              zqns(:,:) 
     1408         IF( ln_mixcpl .AND. ( srcv(jpr_qnsoce)%laction .OR. srcv(jpr_qnsmix)%laction )) THEN  
     1409                                          qns(:,:) = qns(:,:) * xcplmask(:,:,0) + zqns(:,:) * zmsk(:,:)  
     1410         ELSE IF( ll_purecpl ) THEN   ;   qns(:,:) = zqns(:,:) 
    10951411         ENDIF 
    10961412 
     
    11011417         ENDIF 
    11021418         IF( ln_dm2dc .AND. ln_cpl )   zqsr(:,:) = sbc_dcy( zqsr )   ! modify qsr to include the diurnal cycle 
    1103          IF( ln_mixcpl ) THEN   ;   qsr(:,:) = qsr(:,:) * xcplmask(:,:,0) + zqsr(:,:) * zmsk(:,:) 
    1104          ELSE                   ;   qsr(:,:) =                              zqsr(:,:) 
     1419         IF( ln_mixcpl .AND. ( srcv(jpr_qsroce)%laction .OR. srcv(jpr_qsrmix)%laction )) THEN  
     1420                                          qsr(:,:) = qsr(:,:) * xcplmask(:,:,0) + zqsr(:,:) * zmsk(:,:)  
     1421         ELSE IF( ll_purecpl ) THEN   ;   qsr(:,:) = zqsr(:,:) 
    11051422         ENDIF 
    11061423         ! 
     
    11131430      ENDIF 
    11141431      ! 
    1115       CALL wrk_dealloc( jpi,jpj, ztx, zty, zmsk, zemp, zqns, zqsr ) 
     1432      CALL wrk_dealloc( jpi,jpj, ztx, zty, ztx2, zty2, zmsk, zemp, zqns, zqsr ) 
    11161433      ! 
    11171434      IF( nn_timing == 1 )  CALL timing_stop('sbc_cpl_rcv') 
     
    17082025      ! 
    17092026      INTEGER ::   ji, jj, jl   ! dummy loop indices 
     2027      INTEGER ::   ikchoix 
    17102028      INTEGER ::   isec, info   ! local integer 
    17112029      REAL(wp) ::   zumax, zvmax 
     
    17362054            ! 
    17372055            SELECT CASE( sn_snd_temp%cldes) 
    1738             CASE( 'oce only'             )   ;   ztmp1(:,:) =   ztmp1(:,:) + rt0 
     2056            CASE( 'oce only'             )   ;   ztmp1(:,:) = (ztmp1(:,:) + rt0) * tmask(:,:,1) 
    17392057            CASE( 'oce and ice'          )   ;   ztmp1(:,:) =   ztmp1(:,:) + rt0 
    17402058               SELECT CASE( sn_snd_temp%clcat ) 
     
    17652083                  ztmp1(:,:) = ztmp1(:,:) + tn_ice(:,:,jl) * a_i(:,:,jl) 
    17662084               ENDDO 
     2085            CASE( 'none'         )       ! nothing to do 
    17672086            CASE default                     ;   CALL ctl_stop( 'sbc_cpl_snd: wrong definition of sn_snd_temp%cldes' ) 
    17682087            END SELECT 
     
    18892208         !                                                  j+1   j     -----V---F 
    18902209         ! surface velocity always sent from T point                     !       | 
    1891          !                                                        j      |   T   U 
     2210         ! [except for HadGEM3]                                   j      |   T   U 
    18922211         !                                                               |       | 
    18932212         !                                                   j    j-1   -I-------| 
     
    19012220            SELECT CASE( TRIM( sn_snd_crt%cldes ) ) 
    19022221            CASE( 'oce only'             )      ! C-grid ==> T 
    1903                DO jj = 2, jpjm1 
    1904                   DO ji = fs_2, fs_jpim1   ! vector opt. 
    1905                      zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj  ,1) ) 
    1906                      zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji  ,jj-1,1) )  
    1907                   END DO 
    1908                END DO 
     2222               IF ( TRIM( sn_snd_crt%clvgrd ) == 'T' ) THEN    
     2223                  DO jj = 2, jpjm1    
     2224                     DO ji = fs_2, fs_jpim1   ! vector opt.    
     2225                        zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj  ,1) )    
     2226                        zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji,jj-1,1) )     
     2227                     END DO    
     2228                  END DO    
     2229               ELSE    
     2230                  ! Temporarily Changed for UKV    
     2231                  DO jj = 2, jpjm1    
     2232                     DO ji = 2, jpim1    
     2233                        zotx1(ji,jj) = un(ji,jj,1)    
     2234                        zoty1(ji,jj) = vn(ji,jj,1)    
     2235                     END DO    
     2236                  END DO    
     2237               ENDIF  
    19092238            CASE( 'weighted oce and ice' )    
    19102239               SELECT CASE ( cp_ice_msh ) 
     
    19302259                  END DO 
    19312260               CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T 
    1932                   DO jj = 2, jpjm1 
    1933                      DO ji = 2, jpim1   ! NO vector opt. 
    1934                         zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   
    1935                         zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   
    1936                         zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     & 
    1937                            &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj) 
    1938                         zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     & 
    1939                            &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj) 
    1940                      END DO 
    1941                   END DO 
     2261                  IF ( TRIM( sn_snd_crt%clvgrd ) == 'T' ) THEN    
     2262                     DO jj = 2, jpjm1    
     2263                        DO ji = 2, jpim1   ! NO vector opt.    
     2264                           zotx1(ji,jj) = 0.5  * ( un(ji,jj,1) + un(ji-1,jj,1) ) * zfr_l(ji,jj)   &       
     2265                                  &       + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &    
     2266                                  &                + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2267                           zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1) + vn(ji,jj-1,1) ) * zfr_l(ji,jj)   &    
     2268                                  &       + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &    
     2269                                  &                + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2270                        END DO    
     2271                     END DO    
     2272#if defined key_cice    
     2273                  ELSE    
     2274                     ! Temporarily Changed for HadGEM3    
     2275                     DO jj = 2, jpjm1    
     2276                        DO ji = 2, jpim1   ! NO vector opt.    
     2277                           zotx1(ji,jj) = (1.0-fr_iu(ji,jj)) * un(ji,jj,1)             &    
     2278                                &              + fr_iu(ji,jj) * 0.5 * ( u_ice(ji,jj-1) + u_ice(ji,jj) )     
     2279                           zoty1(ji,jj) = (1.0-fr_iv(ji,jj)) * vn(ji,jj,1)             &    
     2280                                &              + fr_iv(ji,jj) * 0.5 * ( v_ice(ji-1,jj) + v_ice(ji,jj) )     
     2281                        END DO    
     2282                     END DO    
     2283#endif    
     2284                  ENDIF  
    19422285               END SELECT 
    19432286               CALL lbc_lnk( zitx1, 'T', -1. )   ;   CALL lbc_lnk( zity1, 'T', -1. ) 
     
    19842327         IF( TRIM( sn_snd_crt%clvor ) == 'eastward-northward' ) THEN             ! Rotation of the components 
    19852328            !                                                                     ! Ocean component 
    1986             CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->e', ztmp1 )       ! 1st component  
    1987             CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->n', ztmp2 )       ! 2nd component  
    1988             zotx1(:,:) = ztmp1(:,:)                                                   ! overwrite the components  
    1989             zoty1(:,:) = ztmp2(:,:) 
    1990             IF( ssnd(jps_ivx1)%laction ) THEN                                     ! Ice component 
    1991                CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 )    ! 1st component  
    1992                CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 )    ! 2nd component  
    1993                zitx1(:,:) = ztmp1(:,:)                                                ! overwrite the components  
    1994                zity1(:,:) = ztmp2(:,:) 
    1995             ENDIF 
     2329             IF ( TRIM( sn_snd_crt%clvgrd ) == 'T' ) THEN    
     2330                CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->e', ztmp1 )       ! 1st component    
     2331                CALL rot_rep( zotx1, zoty1, ssnd(jps_ocx1)%clgrid, 'ij->n', ztmp2 )       ! 2nd component    
     2332                zotx1(:,:) = ztmp1(:,:)                                                   ! overwrite the components    
     2333                zoty1(:,:) = ztmp2(:,:)    
     2334                IF( ssnd(jps_ivx1)%laction ) THEN                  ! Ice component    
     2335                   CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 )    ! 1st component    
     2336                   CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 )    ! 2nd component    
     2337                   zitx1(:,:) = ztmp1(:,:) ! overwrite the components    
     2338                   zity1(:,:) = ztmp2(:,:)    
     2339                ENDIF    
     2340             ELSE    
     2341                ! Temporary code for HadGEM3 - will be removed eventually.    
     2342                ! Only applies when we want uvel on U grid and vvel on V grid    
     2343                ! Rotate U and V onto geographic grid before sending.    
     2344              
     2345                DO jj=2,jpjm1    
     2346                   DO ji=2,jpim1    
     2347                      ztmp1(ji,jj)=0.25*vmask(ji,jj,1)      &    
     2348                           *(zotx1(ji,jj)+zotx1(ji-1,jj)    &    
     2349                           +zotx1(ji,jj+1)+zotx1(ji-1,jj+1))    
     2350                      ztmp2(ji,jj)=0.25*umask(ji,jj,1)      &    
     2351                           *(zoty1(ji,jj)+zoty1(ji+1,jj)    &    
     2352                           +zoty1(ji,jj-1)+zoty1(ji+1,jj-1))    
     2353                   ENDDO    
     2354                ENDDO    
     2355                    
     2356                ! Ensure any N fold and wrap columns are updated    
     2357                CALL lbc_lnk(ztmp1, 'V', -1.0)    
     2358                CALL lbc_lnk(ztmp2, 'U', -1.0)    
     2359                    
     2360                ikchoix = -1    
     2361                CALL repcmo(zotx1,ztmp2,ztmp1,zoty1,zotx1,zoty1,ikchoix)    
     2362            ENDIF    
    19962363         ENDIF 
    19972364         ! 
     
    20192386      ENDIF 
    20202387      ! 
     2388      !                                                      ! ------------------------- !    
     2389      !                                                      !  Surface current to waves !    
     2390      !                                                      ! ------------------------- !    
     2391      IF( ssnd(jps_ocxw)%laction .OR. ssnd(jps_ocyw)%laction ) THEN    
     2392          !        
     2393          !                                                  j+1  j     -----V---F    
     2394          ! surface velocity always sent from T point                    !       |    
     2395          !                                                       j      |   T   U    
     2396          !                                                              |       |    
     2397          !                                                   j   j-1   -I-------|    
     2398          !                                               (for I)        |       |    
     2399          !                                                             i-1  i   i    
     2400          !                                                              i      i+1 (for I)    
     2401          SELECT CASE( TRIM( sn_snd_crtw%cldes ) )    
     2402          CASE( 'oce only'             )      ! C-grid ==> T    
     2403             DO jj = 2, jpjm1    
     2404                DO ji = fs_2, fs_jpim1   ! vector opt.    
     2405                   zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj  ,1) )    
     2406                   zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji , jj-1,1) )     
     2407                END DO    
     2408             END DO    
     2409          CASE( 'weighted oce and ice' )       
     2410             SELECT CASE ( cp_ice_msh )    
     2411             CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T    
     2412                DO jj = 2, jpjm1    
     2413                   DO ji = fs_2, fs_jpim1   ! vector opt.    
     2414                      zotx1(ji,jj) = 0.5 * ( un   (ji,jj,1) + un (ji-1,jj  ,1) ) * zfr_l(ji,jj)      
     2415                      zoty1(ji,jj) = 0.5 * ( vn   (ji,jj,1) + vn (ji  ,jj-1,1) ) * zfr_l(ji,jj)    
     2416                      zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj  ) + u_ice(ji-1,jj    ) ) *  fr_i(ji,jj)    
     2417                      zity1(ji,jj) = 0.5 * ( v_ice(ji,jj  ) + v_ice(ji  ,jj-1  ) ) *  fr_i(ji,jj)    
     2418                   END DO    
     2419                END DO    
     2420             CASE( 'I' )                      ! Ocean on C grid, Ice on I-point (B-grid) ==> T    
     2421                DO jj = 2, jpjm1    
     2422                   DO ji = 2, jpim1   ! NO vector opt.    
     2423                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)      
     2424                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)      
     2425                      zitx1(ji,jj) = 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1)                     &    
     2426                         &                  + u_ice(ji+1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2427                      zity1(ji,jj) = 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1)                     &    
     2428                         &                  + v_ice(ji+1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2429                   END DO    
     2430                END DO    
     2431             CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T    
     2432                DO jj = 2, jpjm1    
     2433                   DO ji = 2, jpim1   ! NO vector opt.    
     2434                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)      
     2435                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)      
     2436                      zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &    
     2437                         &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2438                      zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &    
     2439                         &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2440                   END DO    
     2441                END DO    
     2442             END SELECT    
     2443             CALL lbc_lnk( zitx1, 'T', -1. )   ;   CALL lbc_lnk( zity1, 'T', -1. )    
     2444          CASE( 'mixed oce-ice'        )    
     2445             SELECT CASE ( cp_ice_msh )    
     2446             CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T    
     2447                DO jj = 2, jpjm1    
     2448                   DO ji = fs_2, fs_jpim1   ! vector opt.    
     2449                      zotx1(ji,jj) = 0.5 * ( un   (ji,jj,1) + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &    
     2450                         &         + 0.5 * ( u_ice(ji,jj  ) + u_ice(ji-1,jj    ) ) *  fr_i(ji,jj)    
     2451                      zoty1(ji,jj) = 0.5 * ( vn   (ji,jj,1) + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   &    
     2452                         &         + 0.5 * ( v_ice(ji,jj  ) + v_ice(ji  ,jj-1  ) ) *  fr_i(ji,jj)    
     2453                   END DO    
     2454                END DO    
     2455             CASE( 'I' )                      ! Ocean on C grid, Ice on I-point (B-grid) ==> T    
     2456                DO jj = 2, jpjm1    
     2457                   DO ji = 2, jpim1   ! NO vector opt.    
     2458                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &       
     2459                         &         + 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1)                     &    
     2460                         &                  + u_ice(ji+1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2461                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   &     
     2462                         &         + 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1)                     &    
     2463                         &                  + v_ice(ji+1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2464                   END DO    
     2465                END DO    
     2466             CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T    
     2467                DO jj = 2, jpjm1    
     2468                   DO ji = 2, jpim1   ! NO vector opt.    
     2469                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &       
     2470                         &         + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &    
     2471                         &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2472                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   &     
     2473                         &         + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &    
     2474                         &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)    
     2475                   END DO    
     2476                END DO    
     2477             END SELECT    
     2478          END SELECT    
     2479         CALL lbc_lnk( zotx1, ssnd(jps_ocxw)%clgrid, -1. )   ; CALL lbc_lnk( zoty1, ssnd(jps_ocyw)%clgrid, -1. )    
     2480         !    
     2481         !    
     2482         IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) THEN             ! Rotation of the components    
     2483         !                                                                        ! Ocean component    
     2484            CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->e', ztmp1 )       ! 1st component     
     2485            CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->n', ztmp2 )       ! 2nd component     
     2486            zotx1(:,:) = ztmp1(:,:)                                                   ! overwrite the components     
     2487            zoty1(:,:) = ztmp2(:,:)     
     2488            IF( ssnd(jps_ivx1)%laction ) THEN                                     ! Ice component    
     2489               CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 )    ! 1st component     
     2490               CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 )    ! 2nd component     
     2491               zitx1(:,:) = ztmp1(:,:)                                                ! overwrite the components     
     2492               zity1(:,:) = ztmp2(:,:)    
     2493            ENDIF    
     2494         ENDIF    
     2495         !    
     2496!         ! spherical coordinates to cartesian -> 2 components to 3 components    
     2497!         IF( TRIM( sn_snd_crtw%clvref ) == 'cartesian' ) THEN    
     2498!            ztmp1(:,:) = zotx1(:,:)                     ! ocean currents    
     2499!            ztmp2(:,:) = zoty1(:,:)    
     2500!            CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 )    
     2501!            !    
     2502!            IF( ssnd(jps_ivx1)%laction ) THEN           ! ice velocities    
     2503!               ztmp1(:,:) = zitx1(:,:)    
     2504!               ztmp1(:,:) = zity1(:,:)    
     2505!               CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 )    
     2506!            ENDIF    
     2507!         ENDIF    
     2508         !    
     2509         IF( ssnd(jps_ocxw)%laction )   CALL cpl_snd( jps_ocxw, isec, RESHAPE ( zotx1, (/jpi,jpj,1/) ), info )   ! ocean x current 1st grid    
     2510         IF( ssnd(jps_ocyw)%laction )   CALL cpl_snd( jps_ocyw, isec, RESHAPE ( zoty1, (/jpi,jpj,1/) ), info )   ! ocean y current 1st grid    
     2511         !     
     2512      ENDIF    
     2513      !    
     2514      IF( ssnd(jps_ficet)%laction ) THEN    
     2515         CALL cpl_snd( jps_ficet, isec, RESHAPE ( fr_i, (/jpi,jpj,1/) ), info )    
     2516      END IF    
     2517      !                                                      ! ------------------------- !    
     2518      !                                                      !   Water levels to waves   !    
     2519      !                                                      ! ------------------------- !    
     2520      IF( ssnd(jps_wlev)%laction ) THEN    
     2521         IF( ln_apr_dyn ) THEN     
     2522            IF( kt /= nit000 ) THEN     
     2523               ztmp1(:,:) = sshb(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) )     
     2524            ELSE     
     2525               ztmp1(:,:) = sshb(:,:)     
     2526            ENDIF     
     2527         ELSE     
     2528            ztmp1(:,:) = sshn(:,:)     
     2529         ENDIF     
     2530         CALL cpl_snd( jps_wlev  , isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )    
     2531      END IF 
    20212532      ! 
    20222533      !  Fields sent by OPA to SAS when doing OPA<->SAS coupling 
Note: See TracChangeset for help on using the changeset viewer.