Changeset 7350


Ignore:
Timestamp:
2016-11-28T13:08:46+01:00 (4 years ago)
Author:
emanuelaclementi
Message:

ticket #1805 step 2: Add in changes from the 2015/dev_r5936_INGV1_WAVE branch

Location:
branches/2016/dev_INGV_UKMO_2016/NEMOGCM
Files:
12 edited
2 copied

Legend:

Unmodified
Added
Removed
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/CONFIG/SHARED/namelist_ref

    r5930 r7350  
    250250                           !     =1 global mean of e-p-r set to zero at each time step 
    251251                           !     =2 annual global mean of e-p-r set to zero 
    252    ln_wave = .false.       !  Activate coupling with wave (either Stokes Drift or Drag coefficient, or both)  (T => fill namsbc_wave) 
    253    ln_cdgw = .false.       !  Neutral drag coefficient read from wave model (T => fill namsbc_wave) 
    254    ln_sdw  = .false.       !  Computation of 3D stokes drift                (T => fill namsbc_wave) 
     252   ln_wave = .false.       !  Activate coupling with wave  (T => fill namsbc_wave) 
     253   ln_cdgw = .false.       !  Neutral drag coefficient read from wave model (T => ln_wave=.true. & fill namsbc_wave) 
     254   ln_sdw  = .false.       !  Read 2D Surf Stokes Drift & Computation of 3D stokes drift (T => ln_wave=.true. & fill namsbc_wave)  
     255   ln_tauoc= .false.       !  Activate ocean stress modified by external wave induced stress (T => ln_wave=.true. & fill namsbc_wave) 
     256   ln_stcor= .false.       !  Activate Stokes Coriolis term (T => ln_wave=.true. & ln_sdw=.true. & fill namsbc_wave) 
    255257   nn_lsm  = 0             !  =0 land/sea mask for input fields is not applied (keep empty land/sea mask filename field) , 
    256258                           !  =1:n number of iterations of land/sea mask application for input fields (fill land/sea mask filename field) 
     
    349351   sn_snd_crt    =       'none'                 ,    'no'    , 'spherical' , 'eastward-northward' ,  'T' 
    350352   sn_snd_co2    =       'coupled'              ,    'no'    ,     ''      ,         ''           ,   '' 
     353   sn_snd_crtw   =       'none'                 ,    'no'    ,     ''      ,         ''           , 'U,V' 
     354   sn_snd_ifrac  =       'none'                 ,    'no'    ,     ''      ,         ''           ,   '' 
     355   sn_snd_wlev   =       'coupled'              ,    'no'    ,     ''      ,         ''           ,   '' 
    351356! receive 
    352357   sn_rcv_w10m   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     
    360365   sn_rcv_cal    =       'coupled'              ,    'no'    ,     ''      ,         ''          ,   '' 
    361366   sn_rcv_co2    =       'coupled'              ,    'no'    ,     ''      ,         ''          ,   '' 
     367   sn_rcv_hsig   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     368   sn_rcv_iceflx =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     369   sn_rcv_mslp   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     370   sn_rcv_phioc  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     371   sn_rcv_sdrfx  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     372   sn_rcv_sdrfy  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     373   sn_rcv_wper   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     374   sn_rcv_wnum   =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     375   sn_rcv_wstrf  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
     376   sn_rcv_wdrag  =       'none'                 ,    'no'    ,     ''      ,         ''          ,   '' 
    362377! 
    363378   nn_cplmodel   =     1     !  Maximum number of models to/from which NEMO is potentialy sending/receiving data 
     
    915930   ln_zdfexp   = .false.   !  time-stepping: split-explicit (T) or implicit (F) time stepping 
    916931   nn_zdfexp   =    3            !  number of sub-timestep for ln_zdfexp=T 
     932   ln_zdfqiao  = .false.   !  Enhanced wave vertical mixing Qiao (2010) (T => ln_wave=.true. & ln_sdw=.true. & fill namsbc_wave) 
    917933/ 
    918934!----------------------------------------------------------------------- 
     
    12371253!              !  file name  ! frequency (hours) ! variable     ! time interp. !  clim   ! 'yearly'/ ! weights  ! rotation ! land/sea mask ! 
    12381254!              !             !  (if <0  months)  !   name       !   (logical)  !  (T/F)  ! 'monthly' ! filename ! pairing  ! filename      ! 
    1239    sn_cdg      =  'cdg_wave' ,        1          , 'drag_coeff' ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
     1255   sn_cdg      =  'sdw_wave' ,        1          , 'drag_coeff' ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    12401256   sn_usd      =  'sdw_wave' ,        1          , 'u_sd2d'     ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    12411257   sn_vsd      =  'sdw_wave' ,        1          , 'v_sd2d'     ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    1242    sn_wn       =  'sdw_wave' ,        1          , 'wave_num'   ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
    1243 ! 
    1244    cn_dir_cdg  = './'  !  root directory for the location of drag coefficient files 
     1258   sn_swh      =  'sdw_wave' ,        1          , 'hs'         ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
     1259   sn_wmp      =  'sdw_wave' ,        1          , 'wmp'        ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
     1260   sn_wnum     =  'sdw_wave' ,        1          , 'wave_num'   ,     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
     1261   sn_tauoc    =  'sdw_wave' ,        1          , 'wave_stress',     .true.   , .false. , 'daily'   ,  ''      , ''       , '' 
     1262! 
     1263   cn_dir  = './'  !  root directory for the location of drag coefficient files 
    12451264/ 
    12461265!----------------------------------------------------------------------- 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/SBC/cpl_oasis3.F90

    r5836 r7350  
    6666   INTEGER                    ::   nsnd         ! total number of fields sent  
    6767   INTEGER                    ::   ncplmodel    ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 
    68    INTEGER, PUBLIC, PARAMETER ::   nmaxfld=50   ! Maximum number of coupling fields 
     68   INTEGER, PUBLIC, PARAMETER ::   nmaxfld=55   ! Maximum number of coupling fields 
    6969   INTEGER, PUBLIC, PARAMETER ::   nmaxcat=5    ! Maximum number of coupling fields 
    7070   INTEGER, PUBLIC, PARAMETER ::   nmaxcpl=5    ! Maximum number of coupling fields 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90

    r5836 r7350  
    6565   LOGICAL , PUBLIC ::   ln_cdgw        !: true if neutral drag coefficient from wave model 
    6666   LOGICAL , PUBLIC ::   ln_sdw         !: true if 3d stokes drift from wave model 
     67   LOGICAL , PUBLIC ::   ln_tauoc       !: true if normalized stress from wave is used 
     68   LOGICAL , PUBLIC ::   ln_stcor       !: true if Stokes-Coriolis term is used 
    6769   ! 
    6870   LOGICAL , PUBLIC ::   ln_icebergs    !: Icebergs 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90

    r5583 r7350  
    737737 
    738738      !! Neutral coefficients at 10m: 
    739       IF( ln_cdgw ) THEN      ! wave drag case 
     739      IF( ln_wave .AND. ln_cdgw ) THEN      ! wave drag case 
    740740         cdn_wave(:,:) = cdn_wave(:,:) + rsmall * ( 1._wp - tmask(:,:,1) ) 
    741741         ztmp0   (:,:) = cdn_wave(:,:) 
     
    783783         END IF 
    784784        
    785          IF( ln_cdgw ) THEN      ! surface wave case 
     785         IF( ln_wave .AND. ln_cdgw ) THEN      ! surface wave case 
    786786            sqrt_Cd = vkarmn / ( vkarmn / sqrt_Cd_n10 - zpsi_m_u )  
    787787            Cd      = sqrt_Cd * sqrt_Cd 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r5836 r7350  
    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 
     
    105106   INTEGER, PARAMETER ::   jpr_e3t1st = 41            ! first T level thickness  
    106107   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 
     108   INTEGER, PARAMETER ::   jpr_mslp   = 43            ! mean sea level pressure  
     109   INTEGER, PARAMETER ::   jpr_hsig   = 44            ! Hsig  
     110   INTEGER, PARAMETER ::   jpr_phioc  = 45            ! Wave=>ocean energy flux  
     111   INTEGER, PARAMETER ::   jpr_sdrftx = 46            ! Stokes drift on grid 1  
     112   INTEGER, PARAMETER ::   jpr_sdrfty = 47            ! Stokes drift on grid 2  
     113   INTEGER, PARAMETER ::   jpr_wper   = 48            ! Mean wave period 
     114   INTEGER, PARAMETER ::   jpr_wnum   = 49            ! Mean wavenumber 
     115   INTEGER, PARAMETER ::   jpr_wstrf  = 50            ! Stress fraction adsorbed by waves 
     116   INTEGER, PARAMETER ::   jpr_wdrag  = 51            ! Neutral surface drag coefficient 
     117   INTEGER, PARAMETER ::   jprcv      = 51            ! total number of fields received   
    108118 
    109119   INTEGER, PARAMETER ::   jps_fice   =  1            ! ice fraction sent to the atmosphere 
     
    135145   INTEGER, PARAMETER ::   jps_e3t1st = 27            ! first level depth (vvl) 
    136146   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 
     147   INTEGER, PARAMETER ::   jps_ficet  = 29            ! total ice fraction   
     148   INTEGER, PARAMETER ::   jps_ocxw   = 30            ! currents on grid 1   
     149   INTEGER, PARAMETER ::   jps_ocyw   = 31            ! currents on grid 2 
     150   INTEGER, PARAMETER ::   jps_wlev   = 32            ! water level  
     151   INTEGER, PARAMETER ::   jpsnd      = 32            ! total number of fields sent  
    138152 
    139153   !                                                         !!** namelist namsbc_cpl ** 
     
    149163   ! Received from the atmosphere                     ! 
    150164   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                         
     165   TYPE(FLD_C) ::   sn_rcv_cal, sn_rcv_iceflx, sn_rcv_co2, sn_rcv_mslp                            
     166   ! Send to waves  
     167   TYPE(FLD_C) ::   sn_snd_ifrac, sn_snd_crtw, sn_snd_wlev  
     168   ! Received from waves  
     169   TYPE(FLD_C) ::   sn_rcv_hsig,sn_rcv_phioc,sn_rcv_sdrfx,sn_rcv_sdrfy,sn_rcv_wper,sn_rcv_wnum,sn_rcv_wstrf,sn_rcv_wdrag 
    152170   ! Other namelist parameters                        ! 
    153171   INTEGER     ::   nn_cplmodel            ! Maximum number of models to/from which NEMO is potentialy sending/receiving data 
     
    161179 
    162180   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:)   ::   albedo_oce_mix     ! ocean albedo sent to atmosphere (mix clear/overcast sky) 
     181 
     182   REAL(wp) ::   rpref = 101000._wp   ! reference atmospheric pressure[N/m2]  
     183   REAL(wp) ::   r1_grau              ! = 1.e0 / (grav * rau0)  
    163184 
    164185   INTEGER , ALLOCATABLE, SAVE, DIMENSION(    :) ::   nrcvinfo           ! OASIS info argument 
     
    179200      !!             ***  FUNCTION sbc_cpl_alloc  *** 
    180201      !!---------------------------------------------------------------------- 
    181       INTEGER :: ierr(3) 
     202      INTEGER :: ierr(4) 
    182203      !!---------------------------------------------------------------------- 
    183204      ierr(:) = 0 
     
    190211      ALLOCATE( xcplmask(jpi,jpj,0:nn_cplmodel) , STAT=ierr(3) ) 
    191212      ! 
     213      IF( .NOT. ln_apr_dyn ) ALLOCATE( ssh_ib(jpi,jpj), ssh_ibb(jpi,jpj), apr(jpi, jpj), STAT=ierr(4) )  
     214 
    192215      sbc_cpl_alloc = MAXVAL( ierr ) 
    193216      IF( lk_mpp            )   CALL mpp_sum ( sbc_cpl_alloc ) 
     
    216239      REAL(wp), POINTER, DIMENSION(:,:) ::   zacs, zaos 
    217240      !! 
    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 
     241      NAMELIST/namsbc_cpl/  sn_snd_temp , sn_snd_alb  , sn_snd_thick , sn_snd_crt   , sn_snd_co2,      &  
     242         &                  sn_rcv_w10m, sn_rcv_taumod, sn_rcv_tau   , sn_rcv_dqnsdt, sn_rcv_qsr,      &  
     243         &                  sn_snd_ifrac, sn_snd_crtw , sn_snd_wlev  , sn_rcv_hsig  , sn_rcv_phioc ,   &  
     244         &                  sn_rcv_sdrfx, sn_rcv_sdrfy, sn_rcv_wper  , sn_rcv_wnum  , sn_rcv_wstrf ,   & 
     245         &                  sn_rcv_wdrag, sn_rcv_qns  , sn_rcv_emp   , sn_rcv_rnf   , sn_rcv_cal   ,   & 
     246         &                  sn_rcv_iceflx,sn_rcv_co2  , nn_cplmodel  , ln_usecplmask, sn_rcv_mslp  
    222247      !!--------------------------------------------------------------------- 
    223248      ! 
     
    260285         WRITE(numout,*)'      sea ice heat fluxes             = ', TRIM(sn_rcv_iceflx%cldes), ' (', TRIM(sn_rcv_iceflx%clcat), ')' 
    261286         WRITE(numout,*)'      atm co2                         = ', TRIM(sn_rcv_co2%cldes   ), ' (', TRIM(sn_rcv_co2%clcat   ), ')' 
     287         WRITE(numout,*)'      significant wave heigth         = ', TRIM(sn_rcv_hsig%cldes  ), ' (', TRIM(sn_rcv_hsig%clcat  ), ')'  
     288         WRITE(numout,*)'      wave to oce energy flux         = ', TRIM(sn_rcv_phioc%cldes ), ' (', TRIM(sn_rcv_phioc%clcat ), ')'  
     289         WRITE(numout,*)'      Surface Stokes drift grid u     = ', TRIM(sn_rcv_sdrfx%cldes ), ' (', TRIM(sn_rcv_sdrfx%clcat ), ')'  
     290         WRITE(numout,*)'      Surface Stokes drift grid v     = ', TRIM(sn_rcv_sdrfy%cldes ), ' (', TRIM(sn_rcv_sdrfy%clcat ), ')'  
     291         WRITE(numout,*)'      Mean wave period                = ', TRIM(sn_rcv_wper%cldes  ), ' (', TRIM(sn_rcv_wper%clcat  ), ')'  
     292         WRITE(numout,*)'      Mean wave number                = ', TRIM(sn_rcv_wnum%cldes  ), ' (', TRIM(sn_rcv_wnum%clcat  ), ')'  
     293         WRITE(numout,*)'      Stress frac adsorbed by waves   = ', TRIM(sn_rcv_wstrf%cldes ), ' (', TRIM(sn_rcv_wstrf%clcat ), ')'  
     294         WRITE(numout,*)'      Neutral surf drag coefficient   = ', TRIM(sn_rcv_wdrag%cldes ), ' (', TRIM(sn_rcv_wdrag%clcat ), ')'  
    262295         WRITE(numout,*)'  sent fields (multiple ice categories)' 
    263296         WRITE(numout,*)'      surface temperature             = ', TRIM(sn_snd_temp%cldes  ), ' (', TRIM(sn_snd_temp%clcat  ), ')' 
    264297         WRITE(numout,*)'      albedo                          = ', TRIM(sn_snd_alb%cldes   ), ' (', TRIM(sn_snd_alb%clcat   ), ')' 
    265298         WRITE(numout,*)'      ice/snow thickness              = ', TRIM(sn_snd_thick%cldes ), ' (', TRIM(sn_snd_thick%clcat ), ')' 
     299         WRITE(numout,*)'      total ice fraction              = ', TRIM(sn_snd_ifrac%cldes ), ' (', TRIM(sn_snd_ifrac%clcat ), ')'  
    266300         WRITE(numout,*)'      surface current                 = ', TRIM(sn_snd_crt%cldes   ), ' (', TRIM(sn_snd_crt%clcat   ), ')' 
    267301         WRITE(numout,*)'                      - referential   = ', sn_snd_crt%clvref  
     
    269303         WRITE(numout,*)'                      - mesh          = ', sn_snd_crt%clvgrd 
    270304         WRITE(numout,*)'      oce co2 flux                    = ', TRIM(sn_snd_co2%cldes   ), ' (', TRIM(sn_snd_co2%clcat   ), ')' 
     305         WRITE(numout,*)'      water level                     = ', TRIM(sn_snd_wlev%cldes  ), ' (', TRIM(sn_snd_wlev%clcat  ), ')'  
     306         WRITE(numout,*)'      mean sea level pressure         = ', TRIM(sn_rcv_mslp%cldes  ), ' (', TRIM(sn_rcv_mslp%clcat  ), ')'  
     307         WRITE(numout,*)'      surface current to waves        = ', TRIM(sn_snd_crtw%cldes  ), ' (', TRIM(sn_snd_crtw%clcat  ), ')'  
     308         WRITE(numout,*)'                      - referential   = ', sn_snd_crtw%clvref  
     309         WRITE(numout,*)'                      - orientation   = ', sn_snd_crtw%clvor  
     310         WRITE(numout,*)'                      - mesh          = ', sn_snd_crtw%clvgrd  
    271311         WRITE(numout,*)'  nn_cplmodel                         = ', nn_cplmodel 
    272312         WRITE(numout,*)'  ln_usecplmask                       = ', ln_usecplmask 
     
    307347      !  
    308348      ! Vectors: change of sign at north fold ONLY if on the local grid 
     349      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 
    309350      IF( TRIM( sn_rcv_tau%clvor ) == 'local grid' )   srcv(jpr_otx1:jpr_itz2)%nsgn = -1. 
    310351       
     
    374415         srcv(jpr_ity1)%clgrid = 'V'                  ! i.e. it is always at U- & V-points for i- & j-comp. resp. 
    375416      ENDIF 
     417      ENDIF 
    376418        
    377419      !                                                      ! ------------------------- ! 
     
    470512      !                                                      ! ------------------------- ! 
    471513      srcv(jpr_co2 )%clname = 'O_AtmCO2'   ;   IF( TRIM(sn_rcv_co2%cldes   ) == 'coupled' )    srcv(jpr_co2 )%laction = .TRUE. 
     514 
     515      !                                                      ! ------------------------- !  
     516      !                                                      ! Mean Sea Level Pressure   !  
     517      !                                                      ! ------------------------- !  
     518      srcv(jpr_mslp)%clname = 'O_MSLP'     ;   IF( TRIM(sn_rcv_mslp%cldes  ) == 'coupled' )    srcv(jpr_mslp)%laction = .TRUE.  
     519 
    472520      !                                                      ! ------------------------- ! 
    473521      !                                                      !   topmelt and botmelt     !    
     
    483531         srcv(jpr_topm:jpr_botm)%laction = .TRUE. 
    484532      ENDIF 
     533      !                                                      ! ------------------------- ! 
     534      !                                                      !      Wave breaking        !     
     535      !                                                      ! ------------------------- !  
     536      srcv(jpr_hsig)%clname  = 'O_Hsigwa'    ! significant wave height 
     537      IF( TRIM(sn_rcv_hsig%cldes  ) == 'coupled' )  THEN 
     538         srcv(jpr_hsig)%laction = .TRUE. 
     539         cpl_hsig = .TRUE. 
     540      ENDIF 
     541      srcv(jpr_phioc)%clname = 'O_PhiOce'    ! wave to ocean energy 
     542      IF( TRIM(sn_rcv_phioc%cldes ) == 'coupled' )  THEN 
     543         srcv(jpr_phioc)%laction = .TRUE. 
     544         cpl_phioc = .TRUE. 
     545      ENDIF 
     546      srcv(jpr_sdrftx)%clname = 'O_Sdrfx'    ! Stokes drift in the u direction 
     547      IF( TRIM(sn_rcv_sdrfx%cldes ) == 'coupled' )  THEN 
     548         srcv(jpr_sdrftx)%laction = .TRUE. 
     549         cpl_sdrftx = .TRUE. 
     550      ENDIF 
     551      srcv(jpr_sdrfty)%clname = 'O_Sdrfy'    ! Stokes drift in the v direction 
     552      IF( TRIM(sn_rcv_sdrfy%cldes ) == 'coupled' )  THEN 
     553         srcv(jpr_sdrfty)%laction = .TRUE. 
     554         cpl_sdrfty = .TRUE. 
     555      ENDIF 
     556      srcv(jpr_wper)%clname = 'O_WPer'       ! mean wave period 
     557      IF( TRIM(sn_rcv_wper%cldes  ) == 'coupled' )  THEN 
     558         srcv(jpr_wper)%laction = .TRUE. 
     559         cpl_wper = .TRUE. 
     560      ENDIF 
     561      srcv(jpr_wnum)%clname = 'O_WNum'       ! mean wave number 
     562      IF( TRIM(sn_rcv_wnum%cldes ) == 'coupled' )  THEN 
     563         srcv(jpr_wnum)%laction = .TRUE. 
     564         cpl_wnum = .TRUE. 
     565      ENDIF 
     566      srcv(jpr_wstrf)%clname = 'O_WStrf'     ! stress fraction adsorbed by the wave 
     567      IF( TRIM(sn_rcv_wstrf%cldes ) == 'coupled' )  THEN 
     568         srcv(jpr_wstrf)%laction = .TRUE. 
     569         cpl_wstrf = .TRUE. 
     570      ENDIF 
     571      srcv(jpr_wdrag)%clname = 'O_WDrag'     ! neutral surface drag coefficient 
     572      IF( TRIM(sn_rcv_wdrag%cldes ) == 'coupled' )  THEN 
     573         srcv(jpr_wdrag)%laction = .TRUE. 
     574         cpl_wdrag = .TRUE. 
     575      ENDIF 
     576      !  
    485577      !                                                      ! ------------------------------- ! 
    486578      !                                                      !   OPA-SAS coupling - rcv by opa !    
     
    637729      !                                                      ! ------------------------- ! 
    638730      ssnd(jps_fice)%clname = 'OIceFrc' 
     731      ssnd(jps_ficet)%clname = 'OIceFrcT'  
    639732      ssnd(jps_hice)%clname = 'OIceTck' 
    640733      ssnd(jps_hsnw)%clname = 'OSnwTck' 
     
    645738      ENDIF 
    646739       
     740      IF (TRIM( sn_snd_ifrac%cldes )  == 'coupled') ssnd(jps_ficet)%laction = .TRUE.  
     741 
    647742      SELECT CASE ( TRIM( sn_snd_thick%cldes ) ) 
    648743      CASE( 'none'         )       ! nothing to do 
     
    665760      ssnd(jps_ocy1)%clname = 'O_OCury1'   ;   ssnd(jps_ivy1)%clname = 'O_IVely1' 
    666761      ssnd(jps_ocz1)%clname = 'O_OCurz1'   ;   ssnd(jps_ivz1)%clname = 'O_IVelz1' 
     762      ssnd(jps_ocxw)%clname = 'O_OCurxw'  
     763      ssnd(jps_ocyw)%clname = 'O_OCuryw'  
    667764      ! 
    668765      ssnd(jps_ocx1:jps_ivz1)%nsgn = -1.   ! vectors: change of the sign at the north fold 
     
    685782      END SELECT 
    686783 
     784      ssnd(jps_ocxw:jps_ocyw)%nsgn = -1.   ! vectors: change of the sign at the north fold  
     785         
     786      IF( sn_snd_crtw%clvgrd == 'U,V' ) THEN  
     787         ssnd(jps_ocxw)%clgrid = 'U' ; ssnd(jps_ocyw)%clgrid = 'V'  
     788      ELSE IF( sn_snd_crtw%clvgrd /= 'T' ) THEN  
     789         CALL ctl_stop( 'sn_snd_crtw%clvgrd must be equal to T' )  
     790      ENDIF  
     791      IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) ssnd(jps_ocxw:jps_ocyw)%nsgn = 1.  
     792      SELECT CASE( TRIM( sn_snd_crtw%cldes ) )  
     793         CASE( 'none'                 )   ; ssnd(jps_ocxw:jps_ocyw)%laction = .FALSE.  
     794         CASE( 'oce only'             )   ; ssnd(jps_ocxw:jps_ocyw)%laction = .TRUE.  
     795         CASE( 'weighted oce and ice' )   !   nothing to do  
     796         CASE( 'mixed oce-ice'        )   ; ssnd(jps_ivx1:jps_ivz1)%laction = .FALSE.  
     797         CASE default   ;   CALL ctl_stop( 'sbc_cpl_init: wrong definition of sn_snd_crtw%cldes' )  
     798      END SELECT  
     799 
    687800      !                                                      ! ------------------------- ! 
    688801      !                                                      !          CO2 flux         ! 
    689802      !                                                      ! ------------------------- ! 
    690803      ssnd(jps_co2)%clname = 'O_CO2FLX' ;  IF( TRIM(sn_snd_co2%cldes) == 'coupled' )    ssnd(jps_co2 )%laction = .TRUE. 
     804 
     805      !                                                      ! ------------------------- !  
     806      !                                                      !     Sea surface height    !  
     807      !                                                      ! ------------------------- !  
     808      ssnd(jps_wlev)%clname = 'O_Wlevel' ;  IF( TRIM(sn_snd_wlev%cldes) == 'coupled' )   ssnd(jps_wlev)%laction = .TRUE.  
    691809 
    692810      !                                                      ! ------------------------------- ! 
     
    783901      IF( ln_dm2dc .AND. ln_cpl .AND. ncpl_qsr_freq /= 86400 )   & 
    784902         &   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 
     903      IF( ln_dm2dc .AND. ln_cpl ) ncpl_qsr_freq = 86400 / ncpl_qsr_freq 
    786904 
    787905      CALL wrk_dealloc( jpi,jpj, zacs, zaos ) 
     
    837955      !!                        emp          upward mass flux [evap. - precip. (- runoffs) (- calving)] (ocean only case) 
    838956      !!---------------------------------------------------------------------- 
     957      USE zdf_oce,  ONLY : ln_zdfqiao 
     958 
     959      IMPLICIT NONE 
     960 
    839961      INTEGER, INTENT(in)           ::   kt          ! ocean model time step index 
    840962      INTEGER, INTENT(in)           ::   k_fsbc      ! frequency of sbc (-> ice model) computation  
     
    9921114      IF( srcv(jpr_co2)%laction )   atm_co2(:,:) = frcv(jpr_co2)%z3(:,:,1) 
    9931115#endif 
     1116      !  
     1117      !                                                      ! ========================= !  
     1118      !                                                      ! Mean Sea Level Pressure   !   (taum)  
     1119      !                                                      ! ========================= !  
     1120      !  
     1121      IF( srcv(jpr_mslp)%laction ) THEN                    ! UKMO SHELF effect of atmospheric pressure on SSH  
     1122          IF( kt /= nit000 )   ssh_ibb(:,:) = ssh_ib(:,:)    !* Swap of ssh_ib fields  
     1123 
     1124          r1_grau = 1.e0 / (grav * rau0)               !* constant for optimization  
     1125          ssh_ib(:,:) = - ( frcv(jpr_mslp)%z3(:,:,1) - rpref ) * r1_grau    ! equivalent ssh (inverse barometer)  
     1126          apr   (:,:) =     frcv(jpr_mslp)%z3(:,:,1)                         !atmospheric pressure  
     1127     
     1128          IF( kt == nit000 ) ssh_ibb(:,:) = ssh_ib(:,:)  ! correct this later (read from restart if possible)  
     1129      END IF  
     1130      ! 
     1131      IF( ln_sdw ) THEN  ! Stokes Drift correction activated 
     1132      !                                                      ! ========================= !  
     1133      !                                                      !       Stokes drift u      ! 
     1134      !                                                      ! ========================= !  
     1135         IF( srcv(jpr_sdrftx)%laction ) zusd2dt(:,:) = frcv(jpr_sdrftx)%z3(:,:,1) 
     1136      ! 
     1137      !                                                      ! ========================= !  
     1138      !                                                      !       Stokes drift v      ! 
     1139      !                                                      ! ========================= !  
     1140         IF( srcv(jpr_sdrfty)%laction ) zvsd2dt(:,:) = frcv(jpr_sdrfty)%z3(:,:,1) 
     1141      ! 
     1142      !                                                      ! ========================= !  
     1143      !                                                      !      Wave mean period     ! 
     1144      !                                                      ! ========================= !  
     1145         IF( srcv(jpr_wper)%laction ) wmp(:,:) = frcv(jpr_wper)%z3(:,:,1) 
     1146      ! 
     1147      !                                                      ! ========================= !  
     1148      !                                                      !  Significant wave height  ! 
     1149      !                                                      ! ========================= !  
     1150         IF( srcv(jpr_hsig)%laction ) swh(:,:) = frcv(jpr_hsig)%z3(:,:,1) 
     1151      ! 
     1152      !                                                      ! ========================= !  
     1153      !                                                      !    Vertical mixing Qiao   ! 
     1154      !                                                      ! ========================= !  
     1155         IF( srcv(jpr_wnum)%laction .AND. ln_zdfqiao ) wnum(:,:) = frcv(jpr_wnum)%z3(:,:,1) 
     1156 
     1157         ! Calculate the 3D Stokes drift both in coupled and not fully uncoupled mode 
     1158         IF( srcv(jpr_sdrftx)%laction .OR. srcv(jpr_sdrfty)%laction .OR. srcv(jpr_wper)%laction & 
     1159                                                                    .OR. srcv(jpr_hsig)%laction ) THEN 
     1160            CALL sbc_stokes() 
     1161            IF( ln_zdfqiao .AND. .NOT. srcv(jpr_wnum)%laction ) CALL sbc_qiao() 
     1162         ENDIF 
     1163         IF( ln_zdfqiao .AND. srcv(jpr_wnum)%laction ) CALL sbc_qiao() 
     1164      ENDIF 
     1165      !                                                      ! ========================= !  
     1166      !                                                      ! Stress adsorbed by waves  ! 
     1167      !                                                      ! ========================= !  
     1168      IF( srcv(jpr_wstrf)%laction .AND. ln_tauoc ) tauoc_wave(:,:) = frcv(jpr_wstrf)%z3(:,:,1) 
     1169 
     1170      !                                                      ! ========================= !  
     1171      !                                                      !   Wave drag coefficient   ! 
     1172      !                                                      ! ========================= !  
     1173      IF( srcv(jpr_wdrag)%laction .AND. ln_cdgw ) cdn_wave(:,:) = frcv(jpr_wdrag)%z3(:,:,1) 
    9941174 
    9951175      !  Fields received by SAS when OASIS coupling 
     
    19842164      ENDIF 
    19852165      ! 
     2166      !                                                      ! ------------------------- !  
     2167      !                                                      !  Surface current to waves !  
     2168      !                                                      ! ------------------------- !  
     2169      IF( ssnd(jps_ocxw)%laction .OR. ssnd(jps_ocyw)%laction ) THEN  
     2170          !      
     2171          !                                                  j+1  j     -----V---F  
     2172          ! surface velocity always sent from T point                    !       |  
     2173          !                                                       j      |   T   U  
     2174          !                                                              |       |  
     2175          !                                                   j   j-1   -I-------|  
     2176          !                                               (for I)        |       |  
     2177          !                                                             i-1  i   i  
     2178          !                                                              i      i+1 (for I)  
     2179          SELECT CASE( TRIM( sn_snd_crtw%cldes ) )  
     2180          CASE( 'oce only'             )      ! C-grid ==> T  
     2181             DO jj = 2, jpjm1  
     2182                DO ji = fs_2, fs_jpim1   ! vector opt.  
     2183                   zotx1(ji,jj) = 0.5 * ( un(ji,jj,1) + un(ji-1,jj  ,1) )  
     2184                   zoty1(ji,jj) = 0.5 * ( vn(ji,jj,1) + vn(ji , jj-1,1) )   
     2185                END DO  
     2186             END DO  
     2187          CASE( 'weighted oce and ice' )     
     2188             SELECT CASE ( cp_ice_msh )  
     2189             CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T  
     2190                DO jj = 2, jpjm1  
     2191                   DO ji = fs_2, fs_jpim1   ! vector opt.  
     2192                      zotx1(ji,jj) = 0.5 * ( un   (ji,jj,1) + un   (ji-1,jj  ,1) ) * zfr_l(ji,jj)    
     2193                      zoty1(ji,jj) = 0.5 * ( vn   (ji,jj,1) + vn   (ji  ,jj-1,1) ) * zfr_l(ji,jj)  
     2194                      zitx1(ji,jj) = 0.5 * ( u_ice(ji,jj  ) + u_ice(ji-1,jj    ) ) *  fr_i(ji,jj)  
     2195                      zity1(ji,jj) = 0.5 * ( v_ice(ji,jj  ) + v_ice(ji  ,jj-1  ) ) *  fr_i(ji,jj)  
     2196                   END DO  
     2197                END DO  
     2198             CASE( 'I' )                      ! Ocean on C grid, Ice on I-point (B-grid) ==> T  
     2199                DO jj = 2, jpjm1  
     2200                   DO ji = 2, jpim1   ! NO vector opt.  
     2201                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)    
     2202                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)    
     2203                      zitx1(ji,jj) = 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1)                     &  
     2204                         &                  + u_ice(ji+1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)  
     2205                      zity1(ji,jj) = 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1)                     &  
     2206                         &                  + v_ice(ji+1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)  
     2207                   END DO  
     2208                END DO  
     2209             CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T  
     2210                DO jj = 2, jpjm1  
     2211                   DO ji = 2, jpim1   ! NO vector opt.  
     2212                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)    
     2213                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)    
     2214                      zitx1(ji,jj) = 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &  
     2215                         &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)  
     2216                      zity1(ji,jj) = 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &  
     2217                         &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)  
     2218                   END DO  
     2219                END DO  
     2220             END SELECT  
     2221             CALL lbc_lnk( zitx1, 'T', -1. )   ;   CALL lbc_lnk( zity1, 'T', -1. )  
     2222          CASE( 'mixed oce-ice'        )  
     2223             SELECT CASE ( cp_ice_msh )  
     2224             CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T  
     2225                DO jj = 2, jpjm1  
     2226                   DO ji = fs_2, fs_jpim1   ! vector opt.  
     2227                      zotx1(ji,jj) = 0.5 * ( un   (ji,jj,1) + un   (ji-1,jj  ,1) ) * zfr_l(ji,jj)   &  
     2228                         &         + 0.5 * ( u_ice(ji,jj  ) + u_ice(ji-1,jj    ) ) *  fr_i(ji,jj)  
     2229                      zoty1(ji,jj) = 0.5 * ( vn   (ji,jj,1) + vn   (ji  ,jj-1,1) ) * zfr_l(ji,jj)   &  
     2230                         &         + 0.5 * ( v_ice(ji,jj  ) + v_ice(ji  ,jj-1  ) ) *  fr_i(ji,jj)  
     2231                   END DO  
     2232                END DO  
     2233             CASE( 'I' )                      ! Ocean on C grid, Ice on I-point (B-grid) ==> T  
     2234                DO jj = 2, jpjm1  
     2235                   DO ji = 2, jpim1   ! NO vector opt.  
     2236                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &     
     2237                         &         + 0.25 * ( u_ice(ji+1,jj+1) + u_ice(ji,jj+1)                     &  
     2238                         &                  + u_ice(ji+1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)  
     2239                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   &   
     2240                         &         + 0.25 * ( v_ice(ji+1,jj+1) + v_ice(ji,jj+1)                     &  
     2241                         &                  + v_ice(ji+1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)  
     2242                   END DO  
     2243                END DO  
     2244             CASE( 'F' )                      ! Ocean on C grid, Ice on F-point (B-grid) ==> T  
     2245                DO jj = 2, jpjm1  
     2246                   DO ji = 2, jpim1   ! NO vector opt.  
     2247                      zotx1(ji,jj) = 0.5  * ( un(ji,jj,1)      + un(ji-1,jj  ,1) ) * zfr_l(ji,jj)   &     
     2248                         &         + 0.25 * ( u_ice(ji-1,jj-1) + u_ice(ji,jj-1)                     &  
     2249                         &                  + u_ice(ji-1,jj  ) + u_ice(ji,jj  )  ) *  fr_i(ji,jj)  
     2250                      zoty1(ji,jj) = 0.5  * ( vn(ji,jj,1)      + vn(ji  ,jj-1,1) ) * zfr_l(ji,jj)   &   
     2251                         &         + 0.25 * ( v_ice(ji-1,jj-1) + v_ice(ji,jj-1)                     &  
     2252                         &                  + v_ice(ji-1,jj  ) + v_ice(ji,jj  )  ) *  fr_i(ji,jj)  
     2253                   END DO  
     2254                END DO  
     2255             END SELECT  
     2256          END SELECT  
     2257         CALL lbc_lnk( zotx1, ssnd(jps_ocxw)%clgrid, -1. )   ;   CALL lbc_lnk( zoty1, ssnd(jps_ocyw)%clgrid, -1. )  
     2258         !  
     2259         !  
     2260         IF( TRIM( sn_snd_crtw%clvor ) == 'eastward-northward' ) THEN             ! Rotation of the components  
     2261         !                                                                        ! Ocean component  
     2262            CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->e', ztmp1 )       ! 1st component   
     2263            CALL rot_rep( zotx1, zoty1, ssnd(jps_ocxw)%clgrid, 'ij->n', ztmp2 )       ! 2nd component   
     2264            zotx1(:,:) = ztmp1(:,:)                                                   ! overwrite the components   
     2265            zoty1(:,:) = ztmp2(:,:)   
     2266            IF( ssnd(jps_ivx1)%laction ) THEN                                     ! Ice component  
     2267               CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->e', ztmp1 )    ! 1st component   
     2268               CALL rot_rep( zitx1, zity1, ssnd(jps_ivx1)%clgrid, 'ij->n', ztmp2 )    ! 2nd component   
     2269               zitx1(:,:) = ztmp1(:,:)                                                ! overwrite the components   
     2270               zity1(:,:) = ztmp2(:,:)  
     2271            ENDIF  
     2272         ENDIF  
     2273         !  
     2274!         ! spherical coordinates to cartesian -> 2 components to 3 components  
     2275!         IF( TRIM( sn_snd_crtw%clvref ) == 'cartesian' ) THEN  
     2276!            ztmp1(:,:) = zotx1(:,:)                     ! ocean currents  
     2277!            ztmp2(:,:) = zoty1(:,:)  
     2278!            CALL oce2geo ( ztmp1, ztmp2, 'T', zotx1, zoty1, zotz1 )  
     2279!            !  
     2280!            IF( ssnd(jps_ivx1)%laction ) THEN           ! ice velocities  
     2281!               ztmp1(:,:) = zitx1(:,:)  
     2282!               ztmp1(:,:) = zity1(:,:)  
     2283!               CALL oce2geo ( ztmp1, ztmp2, 'T', zitx1, zity1, zitz1 )  
     2284!            ENDIF  
     2285!         ENDIF  
     2286         !  
     2287         IF( ssnd(jps_ocxw)%laction )   CALL cpl_snd( jps_ocxw, isec, RESHAPE ( zotx1, (/jpi,jpj,1/) ), info )   ! ocean x current 1st grid  
     2288         IF( ssnd(jps_ocyw)%laction )   CALL cpl_snd( jps_ocyw, isec, RESHAPE ( zoty1, (/jpi,jpj,1/) ), info )   ! ocean y current 1st grid  
     2289         !   
     2290      ENDIF  
     2291      !  
     2292      IF( ssnd(jps_ficet)%laction ) THEN  
     2293         CALL cpl_snd( jps_ficet, isec, RESHAPE ( fr_i, (/jpi,jpj,1/) ), info )  
     2294      END IF  
     2295      !                                                      ! ------------------------- !  
     2296      !                                                      !   Water levels to waves   !  
     2297      !                                                      ! ------------------------- !  
     2298      IF( ssnd(jps_wlev)%laction ) THEN  
     2299         IF( ln_apr_dyn ) THEN   
     2300            IF( kt /= nit000 ) THEN   
     2301               ztmp1(:,:) = sshb(:,:) - 0.5 * ( ssh_ib(:,:) + ssh_ibb(:,:) )   
     2302            ELSE   
     2303               ztmp1(:,:) = sshb(:,:)   
     2304            ENDIF   
     2305         ELSE   
     2306            ztmp1(:,:) = sshn(:,:)   
     2307         ENDIF   
     2308         CALL cpl_snd( jps_wlev  , isec, RESHAPE ( ztmp1, (/jpi,jpj,1/) ), info )  
     2309      END IF  
    19862310      ! 
    19872311      !  Fields sent by OPA to SAS when doing OPA<->SAS coupling 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r5836 r7350  
    8888         &             ln_blk_mfs, ln_apr_dyn, nn_ice, nn_ice_embd, ln_dm2dc   , ln_rnf   ,   & 
    8989         &             ln_ssr    , nn_isf    , nn_fwb, ln_cdgw    , ln_wave    , ln_sdw   ,   & 
    90          &             nn_lsm    , nn_limflx , nn_components, ln_cpl 
     90         &             ln_tauoc  , ln_stcor  , nn_lsm, nn_limflx  , nn_components, ln_cpl  
    9191      INTEGER  ::   ios 
    9292      INTEGER  ::   ierr, ierr0, ierr1, ierr2, ierr3, jpm 
     
    131131         WRITE(numout,*) '              ocean-atmosphere coupled formulation       ln_cpl      = ', ln_cpl 
    132132         WRITE(numout,*) '              forced-coupled mixed formulation           ln_mixcpl   = ', ln_mixcpl 
     133         WRITE(numout,*) '              wave physics                               ln_wave     = ', ln_wave 
     134         WRITE(numout,*) '                 Stokes drift corr. to vert. velocity    ln_sdw      = ', ln_sdw 
     135         WRITE(numout,*) '                 wave modified ocean stress              ln_tauoc    = ', ln_tauoc 
     136         WRITE(numout,*) '                 Stokes coriolis term                    ln_stcor    = ', ln_stcor 
     137         WRITE(numout,*) '                 neutral drag coefficient (CORE, MFS)    ln_cdgw     = ', ln_cdgw 
    133138         WRITE(numout,*) '              OASIS coupling (with atm or sas)           lk_oasis    = ', lk_oasis 
    134139         WRITE(numout,*) '              components of your executable            nn_components = ', nn_components 
     
    216221      IF ( ln_wave ) THEN 
    217222      !Activated wave module but neither drag nor stokes drift activated 
    218          IF ( .NOT.(ln_cdgw .OR. ln_sdw) )   THEN 
    219             CALL ctl_warn( 'Ask for wave coupling but nor drag coefficient (ln_cdgw=F) neither stokes drift activated (ln_sdw=F)' ) 
     223         IF ( .NOT.(ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor ) )   THEN 
     224            CALL ctl_warn( 'Ask for wave coupling but ln_cdgw=F, ln_sdw=F, ln_tauoc=F, ln_stcor=F') 
    220225      !drag coefficient read from wave model definable only with mfs bulk formulae and core  
    221226         ELSEIF (ln_cdgw .AND. .NOT.(ln_blk_mfs .OR. ln_blk_core) )       THEN        
    222227             CALL ctl_stop( 'drag coefficient read from wave model definable only with mfs bulk formulae and core') 
     228         ELSEIF (ln_stcor .AND. .NOT. ln_sdw)                             THEN 
     229             CALL ctl_stop( 'Stokes-Coriolis term calculated only if activated Stokes Drift ln_sdw=T') 
    223230         ENDIF 
    224231      ELSE 
    225       IF ( ln_cdgw .OR. ln_sdw  )                                                           &  
     232      IF ( ln_cdgw .OR. ln_sdw .OR. ln_tauoc .OR. ln_stcor )                &  
    226233         &   CALL ctl_stop( 'Not Activated Wave Module (ln_wave=F) but asked coupling ',    & 
    227          &                  'with drag coefficient (ln_cdgw =T) or Stokes drift (ln_sdw=T) ') 
     234         &                  'with drag coefficient (ln_cdgw =T) '  ,                        & 
     235         &                  'or Stokes Drift (ln_sdw=T) ' ,                                 & 
     236         &                  'or ocean stress modification due to waves (ln_tauoc=T) ',      &   
     237         &                  'or Stokes-Coriolis term (ln_stcori=T)'  ) 
    228238      ENDIF  
    229239      !                          ! Choice of the Surface Boudary Condition (set nsbc) 
     
    360370                             CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice )   ! OPA-SAS coupling: OPA receiving fields from SAS 
    361371      END SELECT 
    362  
     372      IF ( ln_wave .AND. ln_tauoc) THEN                                 ! Wave stress subctracted 
     373            utau(:,:) = utau(:,:)*tauoc_wave(:,:) 
     374            vtau(:,:) = vtau(:,:)*tauoc_wave(:,:) 
     375            taum(:,:) = taum(:,:)*tauoc_wave(:,:) 
     376      ! 
     377            SELECT CASE( nsbc ) 
     378            CASE(  0,1,2,3,5,-1 )  ; 
     379                IF(lwp) WRITE(numout,*) 'WARNING: You are subtracting the wave stress to the ocean. & 
     380                        & If not requested select ln_tauoc=.false' 
     381            END SELECT 
     382      ! 
     383      END IF 
    363384      IF( ln_mixcpl )        CALL sbc_cpl_rcv ( kt, nn_fsbc, nn_ice )   ! forced-coupled mixed formulation after forcing 
    364385 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/SBC/sbcwave.F90

    r5860 r7350  
    44   !! Wave module  
    55   !!====================================================================== 
    6    !! History :  3.3  !   2011-09  (Adani M)  Original code: Drag Coefficient  
    7    !!         :  3.4  !   2012-10  (Adani M)                 Stokes Drift  
    8    !!---------------------------------------------------------------------- 
    9  
    10    !!---------------------------------------------------------------------- 
    11    !!   sbc_wave      : read drag coefficient from wave model in netcdf files  
     6   !! History :  3.3  !   2011-09  (M. Adani)  Original code: Drag Coefficient  
     7   !!         :  3.4  !   2012-10  (M. Adani)  Stokes Drift  
     8   !!            3.6  !   2014-09  (E. Clementi,P. Oddo) New Stokes Drift Computation 
     9   !!---------------------------------------------------------------------- 
     10 
     11   !!---------------------------------------------------------------------- 
     12   !!   sbc_wave      : wave data from wave model in netcdf files  
    1213   !!---------------------------------------------------------------------- 
    1314   USE oce            !  
    14    USE sbc_oce        ! Surface boundary condition: ocean fields 
     15   USE sbc_oce       ! Surface boundary condition: ocean fields 
    1516   USE bdy_oce        ! 
    1617   USE domvvl         ! 
    17    ! 
    1818   USE iom            ! I/O manager library 
    1919   USE in_out_manager ! I/O manager 
    2020   USE lib_mpp        ! distribued memory computing library 
    21    USE fldread        ! read input fields 
     21   USE fldread       ! read input fields 
    2222   USE wrk_nemo       ! 
     23   USE phycst         ! physical constants  
    2324 
    2425   IMPLICIT NONE 
    2526   PRIVATE 
    2627 
    27    PUBLIC   sbc_wave    ! routine called in sbc_blk_core or sbc_blk_mfs 
     28   PUBLIC   sbc_stokes, sbc_qiao  ! routines called in sbccpl 
     29   PUBLIC   sbc_wave    ! routine called in sbcmod 
    2830    
    29    INTEGER , PARAMETER ::   jpfld  = 3   ! maximum number of files to read for srokes drift 
    30    INTEGER , PARAMETER ::   jp_usd = 1   ! index of stokes drift  (i-component) (m/s)    at T-point 
    31    INTEGER , PARAMETER ::   jp_vsd = 2   ! index of stokes drift  (j-component) (m/s)    at T-point 
    32    INTEGER , PARAMETER ::   jp_wn  = 3   ! index of wave number                 (1/m)    at T-point 
    33  
    34    TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_cd    ! structure of input fields (file informations, fields read) Drag Coefficient 
    35    TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_sd    ! structure of input fields (file informations, fields read) Stokes Drift 
    36  
    37    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION (:,:)   :: cdn_wave  
    38    REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION (:,:,:) :: usd3d, vsd3d, wsd3d  
    39    REAL(wp),         ALLOCATABLE, DIMENSION (:,:)   :: usd2d, vsd2d, uwavenum, vwavenum  
     31   ! Variables checking if the wave parameters are coupled (if not, they are read from file) 
     32   LOGICAL, PUBLIC     ::   cpl_hsig=.FALSE. 
     33   LOGICAL, PUBLIC     ::   cpl_phioc=.FALSE. 
     34   LOGICAL, PUBLIC     ::   cpl_sdrftx=.FALSE. 
     35   LOGICAL, PUBLIC     ::   cpl_sdrfty=.FALSE. 
     36   LOGICAL, PUBLIC     ::   cpl_wper=.FALSE. 
     37   LOGICAL, PUBLIC     ::   cpl_wnum=.FALSE. 
     38   LOGICAL, PUBLIC     ::   cpl_wstrf=.FALSE. 
     39   LOGICAL, PUBLIC     ::   cpl_wdrag=.FALSE. 
     40 
     41   INTEGER ::   jpfld                ! number of files to read for stokes drift 
     42   INTEGER ::   jp_usd               ! index of stokes drift  (i-component) (m/s)    at T-point 
     43   INTEGER ::   jp_vsd               ! index of stokes drift  (j-component) (m/s)    at T-point 
     44   INTEGER ::   jp_swh               ! index of significant wave hight      (m)      at T-point 
     45   INTEGER ::   jp_wmp               ! index of mean wave period            (s)      at T-point 
     46 
     47   TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_cd    ! structure of input fields (file informations, fields read) Drag Coefficient 
     48   TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_sd    ! structure of input fields (file informations, fields read) Stokes Drift 
     49   TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_wn    ! structure of input fields (file informations, fields read) wave number for Qiao 
     50   TYPE(FLD), ALLOCATABLE, DIMENSION(:)  :: sf_tauoc ! structure of input fields (file informations, fields read) normalized wave stress into the ocean 
     51   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: cdn_wave  
     52   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: swh,wmp, wnum 
     53   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: tauoc_wave 
     54   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: tsd2d 
     55   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:)       :: zusd2dt, zvsd2dt 
     56   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)     :: usd3d, vsd3d, wsd3d  
     57   REAL(wp), PUBLIC, ALLOCATABLE, DIMENSION(:,:,:)     :: usd3dt, vsd3dt 
    4058 
    4159   !! * Substitutions 
     
    4967CONTAINS 
    5068 
     69   SUBROUTINE sbc_stokes( ) 
     70      !!--------------------------------------------------------------------- 
     71      !!                     ***  ROUTINE sbc_stokes  *** 
     72      !! 
     73      !! ** Purpose :   compute the 3d Stokes Drift according to Breivik et al., 
     74      !!                2014 (DOI: 10.1175/JPO-D-14-0020.1) 
     75      !! 
     76      !! ** Method  : - Calculate Stokes transport speed  
     77      !!              - Calculate horizontal divergence  
     78      !!              - Integrate the horizontal divergenze from the bottom  
     79      !! ** action   
     80      !!--------------------------------------------------------------------- 
     81      INTEGER                ::   jj,ji,jk  
     82      REAL(wp)                       ::  ztransp, zfac, zsp0, zk, zus, zvs 
     83      REAL(wp), DIMENSION(:,:,:), POINTER :: ze3hdiv   ! 3D workspace 
     84      !!--------------------------------------------------------------------- 
     85      ! 
     86 
     87      CALL wrk_alloc( jpi,jpj,jpk, ze3hdiv ) 
     88      DO jk = 1, jpk 
     89         DO jj = 1, jpj 
     90            DO ji = 1, jpi 
     91               ! On T grid 
     92               ! Stokes transport speed estimated from Hs and Tmean 
     93               ztransp = 2.0_wp*rpi*swh(ji,jj)**2.0_wp/(16.0_wp*MAX(wmp(ji,jj),0.0000001_wp)) 
     94               ! Stokes surface speed 
     95               zsp0 = SQRT( zusd2dt(ji,jj)**2 + zvsd2dt(ji,jj)**2) 
     96               ! Wavenumber scale 
     97               zk = ABS(zsp0)/MAX(ABS(5.97_wp*ztransp),0.0000001_wp) 
     98               ! Depth attenuation 
     99               zfac = EXP(-2.0_wp*zk*fsdept(ji,jj,jk))/(1.0_wp+8.0_wp*zk*fsdept(ji,jj,jk)) 
     100               ! 
     101               usd3dt(ji,jj,jk) = zfac * zusd2dt(ji,jj) * tmask(ji,jj,jk) 
     102               vsd3dt(ji,jj,jk) = zfac * zvsd2dt(ji,jj) * tmask(ji,jj,jk) 
     103            END DO 
     104         END DO 
     105      END DO  
     106      ! Into the U and V Grid 
     107      DO jk = 1, jpkm1 
     108         DO jj = 1, jpjm1 
     109            DO ji = 1, fs_jpim1 
     110               usd3d(ji,jj,jk) = 0.5 *  umask(ji,jj,jk) *   & 
     111                               &  ( usd3dt(ji,jj,jk) + usd3dt(ji+1,jj,jk) ) 
     112               vsd3d(ji,jj,jk) = 0.5 *  vmask(ji,jj,jk) *   & 
     113                               &  ( vsd3dt(ji,jj,jk) + vsd3dt(ji,jj+1,jk) ) 
     114            END DO 
     115         END DO 
     116      END DO 
     117      ! 
     118      CALL lbc_lnk( usd3d(:,:,:), 'U', -1. ) 
     119      CALL lbc_lnk( vsd3d(:,:,:), 'V', -1. ) 
     120      ! 
     121      DO jk = 1, jpkm1               ! Horizontal divergence 
     122         DO jj = 2, jpj 
     123            DO ji = fs_2, jpi 
     124               ze3hdiv(ji,jj,jk) = (  e2u(ji  ,jj) * usd3d(ji  ,jj,jk)     & 
     125                  &                 - e2u(ji-1,jj) * usd3d(ji-1,jj,jk)     & 
     126                  &                 + e1v(ji,jj  ) * vsd3d(ji,jj  ,jk)     & 
     127                  &                 - e1v(ji,jj-1) * vsd3d(ji,jj-1,jk)   ) * r1_e1e2t(ji,jj) 
     128            END DO 
     129         END DO 
     130      END DO 
     131      ! 
     132      IF( .NOT. AGRIF_Root() ) THEN 
     133         IF( nbondi ==  1 .OR. nbondi == 2 )   ze3hdiv(nlci-1,   :  ,:) = 0._wp      ! east 
     134         IF( nbondi == -1 .OR. nbondi == 2 )   ze3hdiv(  2   ,   :  ,:) = 0._wp      ! west 
     135         IF( nbondj ==  1 .OR. nbondj == 2 )   ze3hdiv(  :   ,nlcj-1,:) = 0._wp      ! north 
     136         IF( nbondj == -1 .OR. nbondj == 2 )   ze3hdiv(  :   ,  2   ,:) = 0._wp      ! south 
     137      ENDIF 
     138      ! 
     139      CALL lbc_lnk( ze3hdiv, 'T', 1. ) 
     140      ! 
     141      DO jk = jpkm1, 1, -1                   ! integrate from the bottom the e3t * hor. divergence 
     142         wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - fse3t_n(:,:,jk) * ze3hdiv(:,:,jk) 
     143      END DO 
     144#if defined key_bdy 
     145      IF( lk_bdy ) THEN 
     146         DO jk = 1, jpkm1 
     147            wsd3d(:,:,jk) = wsd3d(:,:,jk) * bdytmask(:,:) 
     148         END DO 
     149      ENDIF 
     150#endif 
     151      CALL wrk_dealloc( jpi,jpj,jpk, ze3hdiv ) 
     152      ! 
     153   END SUBROUTINE sbc_stokes 
     154 
     155   SUBROUTINE sbc_qiao 
     156      !!--------------------------------------------------------------------- 
     157      !!                     ***  ROUTINE sbc_qiao  *** 
     158      !! 
     159      !! ** Purpose :   Qiao formulation for wave enhanced turbulence 
     160      !!                2010 (DOI: 10.1007/s10236-010-0326)  
     161      !! 
     162      !! ** Method  : -  
     163      !! ** action   
     164      !!--------------------------------------------------------------------- 
     165      INTEGER :: jj, ji 
     166 
     167      ! Calculate the module of the stokes drift on T grid 
     168      !------------------------------------------------- 
     169      DO jj = 1, jpj 
     170         DO ji = 1, jpi 
     171            tsd2d(ji,jj) = SQRT( zusd2dt(ji,jj) * zusd2dt(ji,jj) + zvsd2dt(ji,jj) * zvsd2dt(ji,jj) ) 
     172         END DO 
     173      END DO 
     174      ! 
     175   END SUBROUTINE sbc_qiao 
     176 
    51177   SUBROUTINE sbc_wave( kt ) 
    52178      !!--------------------------------------------------------------------- 
    53       !!                     ***  ROUTINE sbc_apr  *** 
    54       !! 
    55       !! ** Purpose :   read drag coefficient from wave model  in netcdf files. 
     179      !!                     ***  ROUTINE sbc_wave  *** 
     180      !! 
     181      !! ** Purpose :   read wave parameters from wave model  in netcdf files. 
    56182      !! 
    57183      !! ** Method  : - Read namelist namsbc_wave 
    58184      !!              - Read Cd_n10 fields in netcdf files  
    59185      !!              - Read stokes drift 2d in netcdf files  
    60       !!              - Read wave number      in netcdf files  
    61       !!              - Compute 3d stokes drift using monochromatic 
    62       !! ** action  :    
    63       !!--------------------------------------------------------------------- 
    64       INTEGER, INTENT( in  ) ::   kt       ! ocean time step 
     186      !!              - Read wave number in netcdf files  
     187      !!              - Compute 3d stokes drift using Breivik et al.,2014 
     188      !!                formulation 
     189      !! ** action   
     190      !!--------------------------------------------------------------------- 
     191      USE zdf_oce,  ONLY : ln_zdfqiao 
     192 
     193      INTEGER, INTENT( in  ) :: kt       ! ocean time step 
    65194      ! 
    66195      INTEGER                ::   ierror   ! return error code 
    67       INTEGER                ::   ifpr, jj,ji,jk  
    68       INTEGER                ::   ios     ! Local integer output status for namelist read 
    69       TYPE(FLD_N), DIMENSION(jpfld) ::   slf_i     ! array of namelist informations on the fields to read 
     196      INTEGER                ::   ifpr 
     197      INTEGER                ::   ios      ! Local integer output status for namelist read 
     198      ! 
    70199      CHARACTER(len=100)     ::  cn_dir                          ! Root directory for location of drag coefficient files 
    71       TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd, sn_wn   ! informations about the fields to be read 
    72       REAL(wp), DIMENSION(:,:,:), POINTER ::   zusd_t, zvsd_t, ze3hdiv   ! 3D workspace 
    73       !! 
    74       NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_wn 
     200      TYPE(FLD_N), ALLOCATABLE, DIMENSION(:) ::   slf_i     ! array of namelist informations on the fields to read 
     201      TYPE(FLD_N)            ::  sn_cdg, sn_usd, sn_vsd,  & 
     202                             &   sn_swh, sn_wmp, sn_wnum, sn_tauoc      ! informations about the fields to be read 
     203      !! 
     204      NAMELIST/namsbc_wave/  sn_cdg, cn_dir, sn_usd, sn_vsd, sn_swh, sn_wmp, sn_wnum, sn_tauoc 
    75205      !!--------------------------------------------------------------------- 
    76206      ! 
     
    87217         IF(lwm) WRITE ( numond, namsbc_wave ) 
    88218         ! 
    89          IF ( ln_cdgw ) THEN 
    90             ALLOCATE( sf_cd(1), STAT=ierror )           !* allocate and fill sf_wave with sn_cdg 
    91             IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
    92             ! 
    93                                    ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1)   ) 
    94             IF( sn_cdg%ln_tint )   ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 
    95             CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
     219         IF( ln_cdgw ) THEN 
     220            IF( .NOT. cpl_wdrag ) THEN 
     221               ALLOCATE( sf_cd(1), STAT=ierror )           !* allocate and fill sf_wave with sn_cdg 
     222               IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
     223               ! 
     224                                      ALLOCATE( sf_cd(1)%fnow(jpi,jpj,1)   ) 
     225               IF( sn_cdg%ln_tint )   ALLOCATE( sf_cd(1)%fdta(jpi,jpj,1,2) ) 
     226               CALL fld_fill( sf_cd, (/ sn_cdg /), cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
     227            ENDIF 
    96228            ALLOCATE( cdn_wave(jpi,jpj) ) 
    97             cdn_wave(:,:) = 0.0 
    98          ENDIF 
    99          IF ( ln_sdw ) THEN 
    100             slf_i(jp_usd) = sn_usd ; slf_i(jp_vsd) = sn_vsd; slf_i(jp_wn) = sn_wn 
    101             ALLOCATE( sf_sd(3), STAT=ierror )           !* allocate and fill sf_wave with sn_cdg 
    102             IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
    103             ! 
    104             DO ifpr= 1, jpfld 
    105                ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 
    106                IF( slf_i(ifpr)%ln_tint )   ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 
    107             END DO 
    108             CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
    109             ALLOCATE( usd2d(jpi,jpj) , vsd2d(jpi,jpj) , uwavenum(jpi,jpj) , vwavenum(jpi,jpj) ) 
     229         ENDIF 
     230 
     231         IF( ln_tauoc ) THEN 
     232            IF( .NOT. cpl_wstrf ) THEN 
     233               ALLOCATE( sf_tauoc(1), STAT=ierror )           !* allocate and fill sf_wave with sn_tauoc 
     234               IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
     235               ! 
     236                                       ALLOCATE( sf_tauoc(1)%fnow(jpi,jpj,1)   ) 
     237               IF( sn_tauoc%ln_tint )  ALLOCATE( sf_tauoc(1)%fdta(jpi,jpj,1,2) ) 
     238               CALL fld_fill( sf_tauoc, (/ sn_tauoc /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 
     239            ENDIF 
     240            ALLOCATE( tauoc_wave(jpi,jpj) ) 
     241         ENDIF 
     242 
     243         IF( ln_sdw ) THEN 
     244            ! Find out how many fields have to be read from file if not coupled 
     245            jpfld=0 
     246            jp_usd=0; jp_vsd=0; jp_swh=0; jp_wmp=0 
     247            IF( .NOT. cpl_sdrftx ) THEN 
     248               jpfld=jpfld+1 
     249               jp_usd=jpfld 
     250            ENDIF 
     251            IF( .NOT. cpl_sdrfty ) THEN 
     252               jpfld=jpfld+1 
     253               jp_vsd=jpfld 
     254            ENDIF 
     255            IF( .NOT. cpl_hsig ) THEN 
     256               jpfld=jpfld+1 
     257               jp_swh=jpfld 
     258            ENDIF 
     259            IF( .NOT. cpl_wper ) THEN 
     260               jpfld=jpfld+1 
     261               jp_wmp=jpfld 
     262            ENDIF 
     263 
     264            ! Read from file only the non-coupled fields  
     265            IF( jpfld > 0 ) THEN 
     266               ALLOCATE( slf_i(jpfld) ) 
     267               IF( jp_usd > 0 ) slf_i(jp_usd) = sn_usd 
     268               IF( jp_vsd > 0 ) slf_i(jp_vsd) = sn_vsd 
     269               IF( jp_swh > 0 ) slf_i(jp_swh) = sn_swh 
     270               IF( jp_wmp > 0 ) slf_i(jp_wmp) = sn_wmp 
     271               ALLOCATE( sf_sd(jpfld), STAT=ierror )           !* allocate and fill sf_sd with stokes drift 
     272               IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable to allocate sf_wave structure' ) 
     273               ! 
     274               DO ifpr= 1, jpfld 
     275                  ALLOCATE( sf_sd(ifpr)%fnow(jpi,jpj,1) ) 
     276                  IF( slf_i(ifpr)%ln_tint )   ALLOCATE( sf_sd(ifpr)%fdta(jpi,jpj,1,2) ) 
     277               END DO 
     278 
     279               CALL fld_fill( sf_sd, slf_i, cn_dir, 'sbc_wave', 'Wave module ', 'namsbc_wave' ) 
     280            ENDIF 
    110281            ALLOCATE( usd3d(jpi,jpj,jpk),vsd3d(jpi,jpj,jpk),wsd3d(jpi,jpj,jpk) ) 
    111             usd3d(:,:,:) = 0._wp   ;   usd2d(:,:) = 0._wp   ;    uwavenum(:,:) = 0._wp 
    112             vsd3d(:,:,:) = 0._wp   ;   vsd2d(:,:) = 0._wp   ;    vwavenum(:,:) = 0._wp 
     282            ALLOCATE( usd3dt(jpi,jpj,jpk),vsd3dt(jpi,jpj,jpk) ) 
     283            ALLOCATE( swh(jpi,jpj), wmp(jpi,jpj) ) 
     284            ALLOCATE( zusd2dt(jpi,jpj), zvsd2dt(jpi,jpj) ) 
     285            usd3d(:,:,:) = 0._wp 
     286            vsd3d(:,:,:) = 0._wp 
    113287            wsd3d(:,:,:) = 0._wp 
    114          ENDIF 
    115       ENDIF 
    116       ! 
    117       IF ( ln_cdgw ) THEN              !==  Neutral drag coefficient  ==! 
     288            IF( ln_zdfqiao ) THEN     !==  Vertical mixing enhancement using Qiao,2010  ==! 
     289               IF( .NOT. cpl_wnum ) THEN 
     290                  ALLOCATE( sf_wn(1), STAT=ierror )           !* allocate and fill sf_wave with sn_wnum 
     291                  IF( ierror > 0 )   CALL ctl_stop( 'STOP', 'sbc_wave: unable toallocate sf_wave structure' ) 
     292                                         ALLOCATE( sf_wn(1)%fnow(jpi,jpj,1)   ) 
     293                  IF( sn_wnum%ln_tint )  ALLOCATE( sf_wn(1)%fdta(jpi,jpj,1,2) ) 
     294                  CALL fld_fill( sf_wn, (/ sn_wnum /), cn_dir, 'sbc_wave', 'Wave module', 'namsbc_wave' ) 
     295               ENDIF 
     296               ALLOCATE( wnum(jpi,jpj),tsd2d(jpi,jpj) ) 
     297            ENDIF 
     298         ENDIF 
     299      ENDIF 
     300      ! 
     301      IF( ln_cdgw .AND. .NOT. cpl_wdrag ) THEN              !==  Neutral drag coefficient  ==! 
    118302         CALL fld_read( kt, nn_fsbc, sf_cd )      ! read from external forcing 
    119303         cdn_wave(:,:) = sf_cd(1)%fnow(:,:,1) 
    120304      ENDIF 
    121       ! 
    122       IF ( ln_sdw )  THEN              !==  Computation of the 3d Stokes Drift  ==! 
     305 
     306      IF( ln_tauoc .AND. .NOT. cpl_wstrf ) THEN             !==  Wave induced stress  ==! 
     307         CALL fld_read( kt, nn_fsbc, sf_tauoc )      !* read wave norm stress from external forcing 
     308         tauoc_wave(:,:) = sf_tauoc(1)%fnow(:,:,1) 
     309      ENDIF 
     310 
     311      IF( ln_sdw )  THEN                         !==  Computation of the 3d Stokes Drift  ==!  
    123312         ! 
    124          CALL fld_read( kt, nn_fsbc, sf_sd )    !* read drag coefficient from external forcing 
     313         ! Read from file only if the field is not coupled 
     314         IF( jpfld > 0 ) THEN 
     315            CALL fld_read( kt, nn_fsbc, sf_sd )      !* read wave parameters from external forcing 
     316            IF( jp_swh > 0 ) swh(:,:)     = sf_sd(jp_swh)%fnow(:,:,1)   ! significant wave height 
     317            IF( jp_wmp > 0 ) wmp(:,:)     = sf_sd(jp_wmp)%fnow(:,:,1)   ! wave mean period 
     318            IF( jp_usd > 0 ) zusd2dt(:,:) = sf_sd(jp_usd)%fnow(:,:,1)   ! 2D zonal Stokes Drift at T point 
     319            IF( jp_vsd > 0 ) zvsd2dt(:,:) = sf_sd(jp_vsd)%fnow(:,:,1)   ! 2D meridional Stokes Drift at T point 
     320         ENDIF 
    125321         ! 
     322         ! Read also wave number if needed, so that it is available in coupling routines 
     323         IF( ln_zdfqiao .AND. .NOT. cpl_wnum ) THEN 
     324            CALL fld_read( kt, nn_fsbc, sf_wn )      !* read wave parameters from external forcing 
     325            wnum(:,:) = sf_wn(1)%fnow(:,:,1) 
     326         ENDIF 
     327            
     328         !==  Computation of the 3d Stokes Drift according to Breivik et al.,2014 
     329         !(DOI: 10.1175/JPO-D-14-0020.1)==!  
    126330         ! 
    127          CALL wrk_alloc( jpi,jpj,jpk,   zusd_t, zvsd_t, ze3hdiv ) 
    128          !                                      !* distribute it on the vertical 
    129          DO jk = 1, jpkm1 
    130             zusd_t(:,:,jk) = sf_sd(jp_usd)%fnow(:,:,1) * EXP( -2._wp * sf_sd(jp_wn)%fnow(:,:,1) * fsdept_n(:,:,jk) ) 
    131             zvsd_t(:,:,jk) = sf_sd(jp_vsd)%fnow(:,:,1) * EXP( -2._wp * sf_sd(jp_wn)%fnow(:,:,1) * fsdept_n(:,:,jk) ) 
    132          END DO 
    133          !                                      !* interpolate the stokes drift from t-point to u- and v-points 
    134          DO jk = 1, jpkm1 
    135             DO jj = 1, jpjm1 
    136                DO ji = 1, jpim1 
    137                    usd3d(ji,jj,jk) = 0.5_wp * ( zusd_t(ji  ,jj,jk) + zusd_t(ji+1,jj,jk) ) * umask(ji,jj,jk) 
    138                    vsd3d(ji,jj,jk) = 0.5_wp * ( zvsd_t(ji  ,jj,jk) + zvsd_t(ji,jj+1,jk) ) * vmask(ji,jj,jk) 
    139                END DO 
    140             END DO 
    141          END DO 
    142          CALL lbc_lnk( usd3d(:,:,:), 'U', -1. ) 
    143          CALL lbc_lnk( vsd3d(:,:,:), 'V', -1. ) 
    144          ! 
    145          DO jk = 1, jpkm1                       !* e3t * Horizontal divergence  ==! 
    146             DO jj = 2, jpjm1 
    147                DO ji = fs_2, fs_jpim1   ! vector opt. 
    148                   ze3hdiv(ji,jj,jk) = (  e2u(ji  ,jj) * fse3u_n(ji  ,jj,jk) * usd3d(ji  ,jj,jk)     & 
    149                      &                 - e2u(ji-1,jj) * fse3u_n(ji-1,jj,jk) * usd3d(ji-1,jj,jk)     & 
    150                      &                 + e1v(ji,jj  ) * fse3v_n(ji,jj  ,jk) * vsd3d(ji,jj  ,jk)     & 
    151                      &                 - e1v(ji,jj-1) * fse3v_n(ji,jj-1,jk) * vsd3d(ji,jj-1,jk)   ) * r1_e1e2t(ji,jj) 
    152                END DO   
    153             END DO   
    154             IF( .NOT. AGRIF_Root() ) THEN 
    155                IF( nbondi ==  1 .OR. nbondi == 2 )   ze3hdiv(nlci-1,   :  ,jk) = 0._wp      ! east 
    156                IF( nbondi == -1 .OR. nbondi == 2 )   ze3hdiv(  2   ,   :  ,jk) = 0._wp      ! west 
    157                IF( nbondj ==  1 .OR. nbondj == 2 )   ze3hdiv(  :   ,nlcj-1,jk) = 0._wp      ! north 
    158                IF( nbondj == -1 .OR. nbondj == 2 )   ze3hdiv(  :   ,  2   ,jk) = 0._wp      ! south 
    159             ENDIF 
    160          END DO 
    161          CALL lbc_lnk( ze3hdiv, 'T', 1. )  
    162          ! 
    163          DO jk = jpkm1, 1, -1                   !* integrate from the bottom the e3t * hor. divergence 
    164             wsd3d(:,:,jk) = wsd3d(:,:,jk+1) - ze3hdiv(:,:,jk) 
    165          END DO 
    166 #if defined key_bdy 
    167          IF( lk_bdy ) THEN 
    168             DO jk = 1, jpkm1 
    169                wsd3d(:,:,jk) = wsd3d(:,:,jk) * bdytmask(:,:) 
    170             END DO 
    171          ENDIF 
    172 #endif 
    173          CALL wrk_dealloc( jpi,jpj,jpk,   zusd_t, zvsd_t, ze3hdiv ) 
    174          !  
     331         ! Calculate only if no necessary fields are coupled, if not calculate later after coupling 
     332         IF( jpfld == 4 ) THEN 
     333            CALL sbc_stokes() 
     334            IF( ln_zdfqiao .AND. .NOT. cpl_wnum ) THEN 
     335               CALL sbc_qiao() 
     336            ENDIF 
     337         ENDIF 
    175338      ENDIF 
    176339      ! 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/TRA/traadv.F90

    r5930 r7350  
    99   !!            3.7  !  2014-05  (G. Madec)  Add 2nd/4th order cases for CEN and FCT schemes  
    1010   !!             -   !  2014-12  (G. Madec) suppression of cross land advection option 
     11   !!            3.6  !  2015-06  (E. Clementi) Addition of Stokes drift in case of wave coupling 
    1112   !!---------------------------------------------------------------------- 
    1213 
     
    3435   USE wrk_nemo       ! Memory Allocation 
    3536   USE timing         ! Timing 
    36  
    37    USE diaptr          ! Poleward heat transport  
     37   USE sbcwave        ! wave module 
     38   USE sbc_oce        ! surface boundary condition: ocean 
     39   USE diaptr         ! Poleward heat transport  
    3840 
    3941   IMPLICIT NONE 
     
    9597      ! 
    9698      !                                          ! set time step 
     99      zun(:,:,:) = 0.0 
     100      zvn(:,:,:) = 0.0 
     101      zwn(:,:,:) = 0.0 
     102      !     
    97103      IF( neuler == 0 .AND. kt == nit000 ) THEN     ! at nit000 
    98104         r2dtra(:) =  rdttra(:)                          ! = rdtra (restarting with Euler time stepping) 
     
    102108      ! 
    103109      !                                         !==  effective transport  ==! 
    104       DO jk = 1, jpkm1 
    105          zun(:,:,jk) = e2u  (:,:) * fse3u(:,:,jk) * un(:,:,jk)                  ! eulerian transport only 
    106          zvn(:,:,jk) = e1v  (:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
    107          zwn(:,:,jk) = e1e2t(:,:)                 * wn(:,:,jk) 
    108       END DO 
     110      IF( ln_wave .AND. ln_sdw )  THEN 
     111         DO jk = 1, jpkm1 
     112            zun(:,:,jk) = e2u(:,:) * fse3u(:,:,jk) *      & 
     113                        &  ( un(:,:,jk) + usd3d(:,:,jk) )                       ! eulerian transport + Stokes Drift 
     114            zvn(:,:,jk) = e1v(:,:) * fse3v(:,:,jk) *      & 
     115                        &  ( vn(:,:,jk) + vsd3d(:,:,jk) ) 
     116            zwn(:,:,jk) = e1e2t(:,:) *                    & 
     117                        &  ( wn(:,:,jk) + wsd3d(:,:,jk) ) 
     118         END DO 
     119      ELSE 
     120         DO jk = 1, jpkm1 
     121            zun(:,:,jk) = e2u  (:,:) * fse3u(:,:,jk) * un(:,:,jk)               ! eulerian transport only 
     122            zvn(:,:,jk) = e1v  (:,:) * fse3v(:,:,jk) * vn(:,:,jk) 
     123            zwn(:,:,jk) = e1e2t(:,:)                 * wn(:,:,jk) 
     124         END DO 
     125      ENDIF 
    109126      ! 
    110127      IF( ln_vvl_ztilde .OR. ln_vvl_layer ) THEN                                ! add z-tilde and/or vvl corrections 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/ZDF/zdf_oce.F90

    r5836 r7350  
    3535   INTEGER , PUBLIC ::   nn_npc      !: non penetrative convective scheme call  frequency 
    3636   INTEGER , PUBLIC ::   nn_npcp     !: non penetrative convective scheme print frequency 
     37   LOGICAL , PUBLIC ::   ln_zdfqiao  !: Enhanced wave vertical mixing Qiao(2010) formulation flag 
    3738 
    3839 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/ZDF/zdfini.F90

    r5836 r7350  
    5151      INTEGER ::   ioptio, ios       ! local integers 
    5252      !! 
    53       NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp,   & 
    54          &              ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp 
     53      NAMELIST/namzdf/ rn_avm0, rn_avt0, nn_avb, nn_havtb, ln_zdfexp, nn_zdfexp,  & 
     54         &        ln_zdfevd, nn_evdm, rn_avevd, ln_zdfnpc, nn_npc, nn_npcp,       & 
     55         &        ln_zdfqiao 
    5556      !!---------------------------------------------------------------------- 
    5657 
     
    8182         WRITE(numout,*) '      npc call  frequency                 nn_npc    = ', nn_npc 
    8283         WRITE(numout,*) '      npc print frequency                 nn_npcp   = ', nn_npcp 
     84         WRITE(numout,*) '      Qiao formulation flag               ln_zdfqiao=', ln_zdfqiao 
    8385      ENDIF 
    8486 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/step.F90

    r5930 r7350  
    2626   !!            3.6  !  2012-07  (J. Simeon, G. Madec. C. Ethe)  Online coarsening of outputs 
    2727   !!            3.6  !  2014-04  (F. Roquet, G. Madec) New equations of state 
     28   !!            3.6  !  2014-10  (E. Clementi, P. Oddo) Add Qiao vertical mixing in case of waves 
    2829   !!            3.7  !  2014-10  (G. Madec)  LDF simplication  
    2930   !!             -   !  2014-12  (G. Madec) remove KPP scheme 
     
    7576      !!              -8- Outputs and diagnostics 
    7677      !!---------------------------------------------------------------------- 
    77       INTEGER ::   jk      ! dummy loop indice 
     78      INTEGER ::   ji,jj,jk ! dummy loop indice 
    7879      INTEGER ::   indic    ! error indicator if < 0 
    7980      INTEGER ::   kcall    ! optional integer argument (dom_vvl_sf_nxt) 
     
    129130                         CALL zdf_bfr( kstp )         ! bottom friction (if quadratic) 
    130131      !                                               ! Vertical eddy viscosity and diffusivity coefficients 
    131       IF( lk_zdfric  )   CALL zdf_ric( kstp )            ! Richardson number dependent Kz 
    132       IF( lk_zdftke  )   CALL zdf_tke( kstp )            ! TKE closure scheme for Kz 
    133       IF( lk_zdfgls  )   CALL zdf_gls( kstp )            ! GLS closure scheme for Kz 
    134       IF( lk_zdfcst  ) THEN                              ! Constant Kz (reset avt, avm[uv] to the background value) 
     132      IF( lk_zdfric  )   CALL zdf_ric ( kstp )             ! Richardson number dependent Kz 
     133      IF( lk_zdftke  )   CALL zdf_tke ( kstp )             ! TKE closure scheme for Kz 
     134      IF( lk_zdfgls  )   CALL zdf_gls ( kstp )             ! GLS closure scheme for Kz 
     135      IF( ln_zdfqiao )   CALL zdf_qiao( kstp )             ! Qiao vertical mixing  
     136      ! 
     137      IF( lk_zdfcst  ) THEN                                ! Constant Kz (reset avt, avm[uv] to the background value) 
    135138         avt (:,:,:) = rn_avt0 * wmask (:,:,:) 
    136139         avmu(:,:,:) = rn_avm0 * wumask(:,:,:) 
     
    217220                         CALL dyn_adv       ( kstp )  ! advection (vector or flux form) 
    218221                         CALL dyn_vor       ( kstp )  ! vorticity term including Coriolis 
     222      IF( ln_wave .AND. ln_sdw .AND. ln_stcor)           & 
     223               &         CALL dyn_stcor     ( kstp )  ! Stokes-Coriolis forcing 
    219224                         CALL dyn_ldf       ( kstp )  ! lateral mixing 
    220225                         CALL dyn_hpg       ( kstp )  ! horizontal gradient of Hydrostatic pressure 
  • branches/2016/dev_INGV_UKMO_2016/NEMOGCM/NEMO/OPA_SRC/step_oce.F90

    r5930 r7350  
    1919   USE sbcapr           ! surface boundary condition: atmospheric pressure 
    2020   USE sbctide          ! Tide initialisation 
     21   USE sbcwave          ! Wave intialisation 
    2122 
    2223   USE traqsr           ! solar radiation penetration      (tra_qsr routine) 
     
    4142   USE dynzdf           ! vertical diffusion               (dyn_zdf routine) 
    4243   USE dynspg           ! surface pressure gradient        (dyn_spg routine) 
     44   USE dynstcor         ! simp. form of Stokes-Coriolis 
    4345 
    4446   USE dynnxt           ! time-stepping                    (dyn_nxt routine) 
     
    7173   USE zdfric           ! Richardson vertical mixing       (zdf_ric routine) 
    7274   USE zdfmxl           ! Mixed-layer depth                (zdf_mxl routine) 
     75   USE zdfqiao          !Qiao module wave induced mixing   (zdf_qiao routine) 
    7376 
    7477   USE zpshde           ! partial step: hor. derivative     (zps_hde routine) 
Note: See TracChangeset for help on using the changeset viewer.